# Population Based Metaheuristic Algorithm Approach for Analysis of Multi-Item Multi-Period Procurement Lot Sizing Problem

Authors: Prasanna Kumar; Mervin Herbert; Srikanth Rao
Journal: Advances in Operations Research
 (2017)

---

## Abstract

This research study focuses on the optimization of multi-item multi-period procurement lot sizing problem for inventory management. Mathematical model is developed which considers different practical constraints like storage space and budget. The aim is to find optimum order quantities of the product so that total cost of inventory is minimized. The NP-hard mathematical model is solved by adopting a novel ant colony optimization approach. Due to lack of benchmark method specified in the literature to assess the performance of the above approach, another metaheuristic based program of genetic algorithm is also employed to solve the problem. The parameters of genetic algorithm model are calibrated using Taguchi method of experiments. The performance of both algorithms is compared using ANOVA analysis with the real time data collected from a valve manufacturing company. It is verified that two methods have not shown any significant difference as far as objective function value is considered. But genetic algorithm is far better than the ACO method when compared on the basis of CPU execution time.


---

## Body

## 1. Introduction

In the present day competitive markets, it is imperative for the organizations to manage their supply chain as efficiently as possible to sustain their market share and improve profitability. Optimized inventory control is the integral part of supply chain management. Harris [1] had suggested a classical inventory model which was later on extended by several research works to incorporate the different practical constraints under different business scenario. Main extensions which had attracted the attention of many researchers were multi-item and multi-period planning inventory model. The multi-period inventory lot sizing scenario with a single product was introduced by Wagner and Whitin in [2], where a dynamic programming solution algorithm was proposed to obtain feasible solutions to the problem.Another extension of EOQ model was worked out by Das et al. [3] by considering the multi-item inventory model with constant demand under the restrictions on storage area, total average shortage cost, and total average inventory investment cost. In another research work, Nenes et al. [4] proposed a model to manage inventory of a product with irregular demand in multiple periods. A nonlinear goal programming model was developed by Panda et al. [5] for lot sizing for a multi-item inventory problem using penalty functions in a decision-making environment. A combination of inventory control and distribution planning was suggested by D. Kim and C. Kim [6] for a multi-period environment. Another additional practical consideration which was suggested into the classical inventory model was the simulation of all unit discount (AUD) and incremental quantity discount (IQD) policy.On this research front, A. K. Maiti and M. Maiti [7] modeled a multi-item inventory control system of items with AUD and IQD and a combination of these discounts. Sana and Chaudhuri [8] extended this work by considering the delayed payments and their effect on the lot sizing.In this research, the multi-item multi-period inventory model is extended to incorporate varying but deterministic periodic demand rate for optimizing the procurement lot size so that total inventory cost is minimized. The model approaches more realistic business situation by considering the all unit discount, overall storage capacity limit, and budget constraint. The aim of the study is to arrive at the optimal order quantities of all items in different periods such that the total inventory cost including the ordering, holding, and the purchasing costs is minimized. The mixed binary integer model that is developed is overly NP-hard. Therefore, two metaheuristic algorithms of GA and ACO have been proposed for its solution. The two algorithms are used to validate and compare the quality of near-optimum solution and ease of their application.Many of recent research works have seen increased use of Taguchi approach and its hybrid methods with other statistical and metaheuristic techniques for optimization due to their versatility and adaptiveness. Hybrid Taguchi-cuckoo search (HTCS) algorithm has been used extensively in the different fields of study. In this research study, Taguchi optimization technique has been employed to enhance the performance of GA by calibrating the parameters [9, 10].The application of AI technique to multi-item multi-period procurement lot sizing problem with novel GA and ACO is main innovation deployed in this research work. The model usability is enhanced by approaching closer to more practical situations considering the storage area and budget constraints and all unit discount scenario. Also very few applications have been sighted in the research literature where Taguchi experimental design has been deployed for the GA parameters optimization.The rest of paper is organized as follows. In the next section, problem under study is explained in detail. Mathematical formulation of considered multi-item multi-period procurement model is presented. Two metaheuristic algorithms GA and ACO and their novelty application for solving our NP-hard problem are discussed in Section3. This section also contains the discussion on Taguchi method application for GA parameter. Computational results are drawn out in Section 4 and two algorithms are compared on the basis of their performance parameters. Section 5 presents the conclusions and directions for future work.


## 2. Problem Definition

The multi-item multi-period procurement lot sizing problem as applied to scenario under discussion is defined here. It is assumed that a company maintains the inventory of many items for satisfying its customers’ requirement under varying but known demand rates. Normal assumption of constant demand rate is foregone to make the model more practical. The known demand rates may change in different periods within a finite planning horizon havingN periods. The initial inventory of all items is the reserve stock. It is assumed that only one order is placed for a particular item in a given period. Lower limit and upper limit of order quantity are also specified as boundary condition. The ordered quantities of items are received in batch sizes where batch is not permitted to be split. All quantity discount regime is considered. Price discount breakpoints are defined so that if an item is ordered in price break quantity level 1, then no discount will be offered. But if an item is ordered in price break 2 then 5% discount will be offered and if an item is ordered in price break 3 then additional 5% discount will be offered and so on. Moreover, the storage space available to hold the total inventory for each period is constrained. There is a restriction on the maximum budget available for the procurement for each period. The aim is to find out the optimum values of order quantities of all the items for each period such that the total cost of the inventory is minimized and the constraints are satisfied. The above defined problem models many real world inventory control systems. The mathematical formulation of the problem is presented in the next section. However, before doing the mathematical formulation, the necessary notations which are to be used are listed here.s: number of components or items orderedt: number of periods or time horizons over which procurement is madep: number of price breaks for consideration of price discountDx,y: requirement or demand rate of the item x at period yKx: ordering cost of the xth component at the beginning of an intervalOx,y: ordering quantity of uth component in interval vz: number of price discount breakpointsbx,z: zth discount breakpoint of xth component z=1,2,…,w(bu,1=0)S: total storage space availableSx: warehouse space occupied by one unit of the xth componentPx,z: purchasing cost of the xth component at the breakpoint “z”Hx: per unit holding cost of the component xPu: purchasing cost of uth component paid at the start of the intervalOC: total ordering costPC: total purchasing costHC: total holding costIC: total inventory costB: total budget available for the procurementZ1: upper band for OxyEx,y: initial inventory of the component “x” in interval “y”Rx: reserve stock for the item xFu,v,w: a binary decision variable, set equal to one if component c is purchased at price breakpoint w in period v, and zero otherwiseLu,v: a binary decision variable, set equal to one if a purchase of a component u is made in period v, and zero otherwiseM: total warehouse spaceMx: warehouse space for xth component


### 2.1. Mathematical Formulation

The total cost of inventory is sum of the ordering cost, holding cost, and purchase cost. This is represented by the following equation:(1)TC=OC+HC+PC,where TC is the total inventory cost, OC is the ordering cost, HC is holding cost, and PC is the purchase cost.The individual components of the total cost are analyzed as follows.


#### 2.1.1. Ordering Cost

Consider thats components are ordered over t time horizons, and total ordering cost OC is obtained as in the following equation:(2)OC=∑x=1s∑y=1tKxLxy.


#### 2.1.2. Holding Cost

Holding cost can be calculated by multiplying the per unit holding cost with the average inventory during that period as shown in the following equation:(3)holding cost=per unit holding cost×average inventory,Average inventory=12×initial inventory during the start of period+residual inventory at the period end.The inventory of the component x at the beginning of time period y+1 can be related to ordered quantity and final inventory during the period y. The beginning inventory of component x in period y+1 is equal to its initial inventory plus the purchased quantity minus its demand, all in interval y, as related in the following equation:(4)Ex,y+1=Ex,y+Ox,v-Dxy.Using (3) and (4),(5)holding cost=∑x=1s∑y=1tHxExy+Oxy-Dxy2.


#### 2.1.3. Purchase Cost

F
x
y
z is a binary decision variable, set equal to one if component x is purchased at price breakpoint z in period y, and zero otherwise. As a result, the total purchasing cost will be as specified in the following equation:(6)Pur. cost=∑x=1s∑y=1t∑z=1pPxzOxyFxyz.


#### 2.1.4. Total Cost

As a result, the complete mathematical model of the inventory control system is as follows.Minimize(7)Z=∑x=1s∑y=1tKxLxy+∑x=1s∑y=1tHxExy+Oxy-Dxy2+∑x=1s∑y=1t∑z=1pPxyOxyFxyz,subject to the constraints:forx=1,2,…,s, y=1,2,…,t, and w=1,2,…,p,(8)1Lx,v=0or1(Boolean value),Fxyz=0or1(Boolean value),(9)2∑z=1pFXYZ=1.Constraint (9) ensures that the quantity should be bought at only one price break.(10)3Ex,y+1=Ex,y+Ox,v-Dxy.Constraint (10) is the mathematical representation of following equation:Inventory brought forward to next period = inventory brought forward to this period + quantity ordered in this period − demand (or consumption) in this period.(11)4Exy≥Kx.Constraint (11) ensures that the inventory brought forward should be greater than or equal to the reserve stock.(12)5∑x=1sPxOxy≤Byfory=1,2,…,p.Equation (12) ensures that budget constraint should be satisfied for each period. (13)6∑x=1smxExy+Oxy≤Mfory=1,2,…,p.Equation (13) ensures that warehouse area constraint should be satisfied for each period.


## 2.1. Mathematical Formulation

The total cost of inventory is sum of the ordering cost, holding cost, and purchase cost. This is represented by the following equation:(1)TC=OC+HC+PC,where TC is the total inventory cost, OC is the ordering cost, HC is holding cost, and PC is the purchase cost.The individual components of the total cost are analyzed as follows.


### 2.1.1. Ordering Cost

Consider thats components are ordered over t time horizons, and total ordering cost OC is obtained as in the following equation:(2)OC=∑x=1s∑y=1tKxLxy.


### 2.1.2. Holding Cost

Holding cost can be calculated by multiplying the per unit holding cost with the average inventory during that period as shown in the following equation:(3)holding cost=per unit holding cost×average inventory,Average inventory=12×initial inventory during the start of period+residual inventory at the period end.The inventory of the component x at the beginning of time period y+1 can be related to ordered quantity and final inventory during the period y. The beginning inventory of component x in period y+1 is equal to its initial inventory plus the purchased quantity minus its demand, all in interval y, as related in the following equation:(4)Ex,y+1=Ex,y+Ox,v-Dxy.Using (3) and (4),(5)holding cost=∑x=1s∑y=1tHxExy+Oxy-Dxy2.


### 2.1.3. Purchase Cost

F
x
y
z is a binary decision variable, set equal to one if component x is purchased at price breakpoint z in period y, and zero otherwise. As a result, the total purchasing cost will be as specified in the following equation:(6)Pur. cost=∑x=1s∑y=1t∑z=1pPxzOxyFxyz.


### 2.1.4. Total Cost

As a result, the complete mathematical model of the inventory control system is as follows.Minimize(7)Z=∑x=1s∑y=1tKxLxy+∑x=1s∑y=1tHxExy+Oxy-Dxy2+∑x=1s∑y=1t∑z=1pPxyOxyFxyz,subject to the constraints:forx=1,2,…,s, y=1,2,…,t, and w=1,2,…,p,(8)1Lx,v=0or1(Boolean value),Fxyz=0or1(Boolean value),(9)2∑z=1pFXYZ=1.Constraint (9) ensures that the quantity should be bought at only one price break.(10)3Ex,y+1=Ex,y+Ox,v-Dxy.Constraint (10) is the mathematical representation of following equation:Inventory brought forward to next period = inventory brought forward to this period + quantity ordered in this period − demand (or consumption) in this period.(11)4Exy≥Kx.Constraint (11) ensures that the inventory brought forward should be greater than or equal to the reserve stock.(12)5∑x=1sPxOxy≤Byfory=1,2,…,p.Equation (12) ensures that budget constraint should be satisfied for each period. (13)6∑x=1smxExy+Oxy≤Mfory=1,2,…,p.Equation (13) ensures that warehouse area constraint should be satisfied for each period.


## 2.1.1. Ordering Cost

Consider thats components are ordered over t time horizons, and total ordering cost OC is obtained as in the following equation:(2)OC=∑x=1s∑y=1tKxLxy.


## 2.1.2. Holding Cost

Holding cost can be calculated by multiplying the per unit holding cost with the average inventory during that period as shown in the following equation:(3)holding cost=per unit holding cost×average inventory,Average inventory=12×initial inventory during the start of period+residual inventory at the period end.The inventory of the component x at the beginning of time period y+1 can be related to ordered quantity and final inventory during the period y. The beginning inventory of component x in period y+1 is equal to its initial inventory plus the purchased quantity minus its demand, all in interval y, as related in the following equation:(4)Ex,y+1=Ex,y+Ox,v-Dxy.Using (3) and (4),(5)holding cost=∑x=1s∑y=1tHxExy+Oxy-Dxy2.


## 2.1.3. Purchase Cost

F
x
y
z is a binary decision variable, set equal to one if component x is purchased at price breakpoint z in period y, and zero otherwise. As a result, the total purchasing cost will be as specified in the following equation:(6)Pur. cost=∑x=1s∑y=1t∑z=1pPxzOxyFxyz.


## 2.1.4. Total Cost

As a result, the complete mathematical model of the inventory control system is as follows.Minimize(7)Z=∑x=1s∑y=1tKxLxy+∑x=1s∑y=1tHxExy+Oxy-Dxy2+∑x=1s∑y=1t∑z=1pPxyOxyFxyz,subject to the constraints:forx=1,2,…,s, y=1,2,…,t, and w=1,2,…,p,(8)1Lx,v=0or1(Boolean value),Fxyz=0or1(Boolean value),(9)2∑z=1pFXYZ=1.Constraint (9) ensures that the quantity should be bought at only one price break.(10)3Ex,y+1=Ex,y+Ox,v-Dxy.Constraint (10) is the mathematical representation of following equation:Inventory brought forward to next period = inventory brought forward to this period + quantity ordered in this period − demand (or consumption) in this period.(11)4Exy≥Kx.Constraint (11) ensures that the inventory brought forward should be greater than or equal to the reserve stock.(12)5∑x=1sPxOxy≤Byfory=1,2,…,p.Equation (12) ensures that budget constraint should be satisfied for each period. (13)6∑x=1smxExy+Oxy≤Mfory=1,2,…,p.Equation (13) ensures that warehouse area constraint should be satisfied for each period.


## 3. Solution Methodology

Since the mathematical model proposed is NP-hard, two metaheuristic approaches ant colony optimization (ACO) and genetic algorithm (GA) have been proposed to solve the problem, with the objective of verifying the solutions obtained and comparing the performances.


### 3.1. Ant Colony Optimization

The phenomenon of pheromone communication has been the drive behind ant colony optimization (ACO) metaheuristics algorithms. Artificial ants simulate the real life ants which communicate their experience while optimizing their search for food in nature through trails of pheromone.Candidate solutions to the optimization problem are constructed by individual ants by interactively adding solution components to initialized empty solution. A complete solution is generated by the ants using the two components: pheromone information which is the accumulated experience and heuristic information which is problem specific data.Which ants are allowed to modify the pheromone information and how they modify is governed by the update strategy. Usually better solution components will receive higher amount of pheromone and will have higher probability of being used by other ants in the subsequent iterations of algorithm.


#### 3.1.1. Ant Cycle Model

Ant cycle model is adopted in this research work. In this model, the trails are globally updated during each cycle by all ants. The amount of pheromone deposited by each ant is a function of the solution quality. As per the flowchart of proposed ant colony optimization algorithm, first the ant based initial solution construction method is executed. Next values of pheromone are set based on ant based initial solution construction method. Thus the values of pheromone are set accurately after the execution of ant based initial solution construction method. Followed by this, ant based solution construction method will be executed in each cycle which will utilize the values of pheromone in order to construct better solutions [11].This research work adopts ant cycle model where the trails are uniformly updated during each cycle by all ants. The solution quality is determined by the amount of pheromone deposited by each ant. The initial solution construction method is executed. Then the ant based solution construction method is executed which will employ the values of pheromone to construct better solutions. Each ant represents a solution to the optimization problem. In our case of application to inventory management, it represents the order quantity of different items in different periods so that total inventory cost is reduced.The novel concept is implemented for ant based solution construction in each cycle as follows. First each ant will construct the solution according to the probabilistic rule. After constructing the solution according to the probabilistic rule, ant 1 will discard its order quantity values of all periods for item 1 and item 3 and the values of order quantity of all periods for item 1 and item 3 of ant 2 will be copied into the solution generated by ant 1. Similarly, ant 2 will discard its order quantity values of all periods for item 1 and item 3 and the values of order quantity of all periods for item 1 and item 3 of ant 3 will be copied into the solution generated by ant 2. This discard-copy step is repeated for all ants. Last ant will discard its order quantity values of all periods for item 1 and item 3 and the values of order quantity of all periods for item 1 and item 3 of ant 1 will be copied into the solution generated by last ant. This novel concept represents the direct communication between ants in each cycle in order to further improve the solution constructed by individual ant separately.In this research work, ant based solution is constructed by using problem specific knowledge so that it does not violate constraint (3) and constraint (4), specified under mathematical modeling under Section2. Once the value of order quantity is generated for all components for all periods, it will be verified that whether the generated values of order quantity will meet the budget constraint and warehouse area constraint or not. If either of the constraints is not met, then the generated values of order quantity (solution) will be discarded. If both the budget constraint and the warehouse area constraint are met, then the solution is valid and it is stored for further processing.The ACO algorithm has been implemented on JAVA platform and the program can be run with Net Beans IDE.


#### 3.1.2. Determination of ACO Parameter

To determine the optimum value for number of cycles, which is an important ACO algorithm parameter, ACO is run several times for one of the data sets, changing the number of cycles each time and recording the objective value function and CPU processing time. ACO run results are tabulated in Table1. From Figure 1 which shows the variation of objective function value with number of cycles in the ACO model, it is clear that objective function reaches optimum value at 600 cycles. After that increasing the number of cycles will only increase the CPU time without improving the objective function value.Table 1
Variation of objective function value and time of execution with number of cycles.
Number of cycles
Objective function value (Rs)
Time of execution (Sec)
50
8382155
5
100
8255916
9
200
8101560
20
300
8094680
32
500
8088430
41
600
8082500
54
700
8082435
62
800
8082440
75
1000
8082440
90Figure 1
Variation of objective function value with number of cycles: ACO model.


### 3.2. Genetic Algorithm

Genetic algorithm is a population based metaheuristic which can be utilized to solve complex optimization problems. John Holland invented genetic algorithm which mimics the principle of natural genetics in order to solve optimization problems. Genetic algorithm uses knowledge from previous generations to direct the future search. A number of chromosomes rather than a single chromosome can be used to carry out the search process. The optimum or near-optimum solution can be efficiently obtained using genetic algorithm [12].The GA implementation comprises the following important stages [13]:(i)
Setting GA parameters including the crossover probability (Pc), the mutation probability (Pm), population size (PS), and number of generations (NG).(ii)
Initializing the population with the size of PS randomly.(iii)
Evaluating the objective function.(iv)
Selecting individual for mating pool by roulette wheel selection method and using elitisms.(v)
Applying the crossover operation for each pair of chromosomes based on Pc.(vi)
Applying mutation operation for each chromosome based on Pm.(vii)
Replacing the current population by the resulting new population.(viii)
Continuing the cycle till stopping criterion is met, like the maximum number of cycles.The GA has been implemented on JAVA platform and the program can be run with Net Beans IDE. The following components of genetic algorithm have been implemented to solve the problem.


#### 3.2.1. Chromosome Representation

The chromosomes represent the strings of the order quantities of the items in each period in this research work. For example, the chromosome representation with four items and three periods is shown in (14). In (14), the number of items and periods are represented by rows and columns, respectively.(14)CR=Q11Q12Q13Q21Q22Q23Q31Q32Q33Q41Q42Q43.


#### 3.2.2. Initial Population Generation

Initial population of size PS is determined using Taguchi’s design which is explained in Section3.3.


#### 3.2.3. Calculate Fitness of Each Chromosome

The process of calculating the fitness of each chromosome consists of the following two steps.(1)
Calculate the objective function of each chromosome.(2)
Calculate the fitness of each chromosome based on its objective function value using the following equation:(15)fx=11+Ox.In (15), fx is the fitness function of each chromosome and Ox is the objective function of each chromosome.


#### 3.2.4. Selection Operator

The roulette wheel selection is used to select the chromosomes. In roulette wheel selection, chromosomes are selected into the mating pool according to their raw fitness. The chromosomes having high fitness have more chances to be selected. The probability ofith-selected chromosome is given by the following equation:(16)Pi=fi∑j=1nfj.In (16), f(i) and f(j) are the fitness of chromosomes i and j, respectively [12].


#### 3.2.5. Crossover Operator

Single point crossover is used in this research work. The following equation (17) demonstrates the application of single point crossover on the two selected parents. Crossover rate is also optimized using Taguchi’s design in order to protect some of good chromosomes that are already present in the mating pool.(17)Parent1=Q11Q12Q13Q21Q22Q23Q31Q32Q33Q41Q42Q43,Offspring1=Q11Q12O13Q21Q22O23Q31Q32O33Q41Q42O43»Parent2=O11O12O13O21O22O23O31O32O33O41O42O43,Offspring2=O11O12Q13O21O22Q23O31O32Q33O41O42Q43.


#### 3.2.6. Mutation Operator

A novel mutation operation, namely, “subtraction-addition” mutation, is used in this research work as follows.Step 1.
Any one period is selected randomly for each component.Step 2.
Mutation value is randomly generated.Step 3.
If selected period ist then mutation value is subtracted from period t and added to period T. If selected period t is the last period then T is the first period. If selected period t is not the last period then T=t+1 period.This novel mutation operation never violates constraint (3) of the total cost optimizing function which is stated in Section2.1.4. Again the mutation value is generated in a way that this novel mutation operation never violates constraint (4) too. Moreover, mutation rate is also optimized using Taguchi design.


#### 3.2.7. Replacement Strategy

The chromosomes in old population are replaced by newly generated chromosomes which meet both the budget constraint (5) and warehouse area constraint (6) in each generation. The newly generated chromosomes which do not meet either constraint (5) or constraint (6) will be discarded.


#### 3.2.8. Stopping Criteria

Genetic algorithm stops when it has completed maximum number of generations. Maximum number of generations can be increased or decreased based on the input size of the problem. This is another parameter of GA which is optimized using Taguchi experiment design.


### 3.3. Parameter Calibration Using Taguchi Design of Experiments

Selection of optimum parameters apparently has significant effect on the effectiveness of the metaheuristic algorithms. When the number of factors affecting the performance of genetic algorithm is more, the full factorial experiment and exhaustive approach to investigate the effect of different parameters becomes increasingly complicated and impractical.In this work, the Taguchi method is used to tune the parameters of metaheuristic algorithm of multiobjective optimization, genetic algorithm. The parameters of a metaheuristic that are needed to be tuned act like controllable factors in the design of experiments (DOE). The aim is to find an optimal combination of the parameters such that the response (the fitness function) is optimized. The Taguchi method is a special case of the fractional factorial design in which some special orthogonal arrays are used [14]. Orthogonal arrays under Taguchi method would help to study a large number of decision variables with a limited number of experiments. Decision variables are divided into controllable and noise factors. Noise factors can not be controlled directly. It is also impractical and most of the time impossible to eliminate the noise factors [15]. Taguchi experimental design will help to reduce the effect of noise factors.Taguchi adopted the concept of signal to noise ratio to reduce the effect of noise factors in the experiment. The desired value or mean response value is represented by signal. The undesirable value or standard deviation is denoted by noise. The variation present in the response variable or the component of noise factor is represented byS/N ratio. Objective functions are classified into 3 types for design of experiment applications by Taguchi. They are “smaller the better,” the “larger the better,” and “the nominal is best.” Since almost all objective functions in inventory control systems are grouped in the “smaller the better” type, its corresponding S/N ratio is given by the following equation:(18)SN=-10log10⁡k2,where k is the objective function value.


#### 3.3.1. Taguchi Method Implementation

Taguchi method implementation consists of five basic steps [14].(i)
Parameters with significant effects on the response are determined.(ii)
The parameter value is determined by the trial and error procedure so that good fitness value is obtained while implementing the experiments.(iii)
Taking into consideration the available degree of freedom (DOF), a suitable orthogonal array is selected which would specify the number of experiments to minimize the experimentation time and at the same time would be giving full weightage for all the combination of influencing factors.(iv)
Experiments are conducted based on obtained design.(v)
The results are recorded. TheS/N approach will give the basis for evaluation and analysis of the results.Setting GA parameters including the crossover probability (Pc), the mutation probability (Pm), population size (PS), and number of generations (NG) is very important in determining the efficiency of the metaheuristic algorithms like GA. Taguchi design which is based on fractional factorial experiments is used to set the GA parameters (see Table2). Minitab 15 is used to employ the Taguchi method. For 3 levels of 4 factors, L9 orthogonal array is suggested which will facilitate representative Fractional Factorial Experimentation. L9 orthogonal array lists the different combinations of factors at different levels at which the response value of experiments has to be determined. 3 items, 3 periods, and 3 price breaks data are selected for experimentation and response in terms of minimum total cost is tabulated as shown in Table 3. For different combinations of the factor levels, each example is solved three times and the mean response was used in the analysis. Figure 2 shows the main effect plot of S/N ratio for different parameter levels of the proposed algorithms. Based on “smaller the better” definition for the signal to noise ratio, optimum values of the GA parameters are as shown in the Table 4.Table 2
GA parameters and levels for Taguchi design.
Levels
GA parameters
Crossover probability Pc
Mutation probability Pm
Population size PS
Number of generations NG
1
0.7
0.15
4000
500
2
0.8
0.20
5000
600
3
0.9
0.30
6000
700Table 3
Experimental response for Taguchi design.
Crossover  probability Pc
Mutation probability Pm
Population size PS
Number of generations NG
Response value obj. function
S
/
N ratio
0.7
0.15
4000
500
8070423
-
138.138
0.7
0.20
5000
600
8071445
-
138.139
0.7
0.30
6000
700
8070334
-
138.138
0.8
0.15
5000
700
8070656
-
138.138
0.8
0.20
6000
500
8072345
-
138.140
0.8
0.30
4000
600
8074558
-
138.142
0.9
0.15
6000
600
8066732
-
138.134
0.9
0.20
4000
700
8065550
-
138.133
0.9
0.30
5000
500
8058597
-
138.125Table 4
Optimum GA parameters.
Optimum GA parameters
Crossover probability Pc
Mutation probability Pm
Population size PS
Number of generations NG
0.9
0.30
5000
500Figure 2
Main effect plot ofS/N ratios: Taguchi design.Minitab 15 is used to employ the Taguchi method. Under the menu options of Minitab, Stat-DOE-Taguchi Design-Create Design is selected. For 3 levels of 4 factors, L9 orthogonal array is suggested which will facilitate representative Fractional Factorial Experimentation. L9 orthogonal array lists the different combinations of factors at different levels at which the response values of experiments have to be determined. 3 items, 3 periods, and 3 price breaks data are selected for experimentation.For different combinations of the factor levels, each example is solved three times and the mean response was used in the analysis. Figure2 shows the main effect plot of S/N ratio and Figure 3 shows the main effect plot of means for different parameter levels of the proposed algorithms.Figure 3
Main effect plot of means: Taguchi design.TheS/N ratio indicates the amount of variation present in the response variable, and the aim is to maximize it. Pc value of 0.9, Pm value of 0.3, PS value of 5000, and NG value of 500 yield the maximum value of S/N ratio as can be seen from Figure 2. This can also be verified from Figure 3 where the above values give the best or the lowest value of objective function. Optimal parameter values of the algorithms are shown in Table 4.


## 3.1. Ant Colony Optimization

The phenomenon of pheromone communication has been the drive behind ant colony optimization (ACO) metaheuristics algorithms. Artificial ants simulate the real life ants which communicate their experience while optimizing their search for food in nature through trails of pheromone.Candidate solutions to the optimization problem are constructed by individual ants by interactively adding solution components to initialized empty solution. A complete solution is generated by the ants using the two components: pheromone information which is the accumulated experience and heuristic information which is problem specific data.Which ants are allowed to modify the pheromone information and how they modify is governed by the update strategy. Usually better solution components will receive higher amount of pheromone and will have higher probability of being used by other ants in the subsequent iterations of algorithm.


### 3.1.1. Ant Cycle Model

Ant cycle model is adopted in this research work. In this model, the trails are globally updated during each cycle by all ants. The amount of pheromone deposited by each ant is a function of the solution quality. As per the flowchart of proposed ant colony optimization algorithm, first the ant based initial solution construction method is executed. Next values of pheromone are set based on ant based initial solution construction method. Thus the values of pheromone are set accurately after the execution of ant based initial solution construction method. Followed by this, ant based solution construction method will be executed in each cycle which will utilize the values of pheromone in order to construct better solutions [11].This research work adopts ant cycle model where the trails are uniformly updated during each cycle by all ants. The solution quality is determined by the amount of pheromone deposited by each ant. The initial solution construction method is executed. Then the ant based solution construction method is executed which will employ the values of pheromone to construct better solutions. Each ant represents a solution to the optimization problem. In our case of application to inventory management, it represents the order quantity of different items in different periods so that total inventory cost is reduced.The novel concept is implemented for ant based solution construction in each cycle as follows. First each ant will construct the solution according to the probabilistic rule. After constructing the solution according to the probabilistic rule, ant 1 will discard its order quantity values of all periods for item 1 and item 3 and the values of order quantity of all periods for item 1 and item 3 of ant 2 will be copied into the solution generated by ant 1. Similarly, ant 2 will discard its order quantity values of all periods for item 1 and item 3 and the values of order quantity of all periods for item 1 and item 3 of ant 3 will be copied into the solution generated by ant 2. This discard-copy step is repeated for all ants. Last ant will discard its order quantity values of all periods for item 1 and item 3 and the values of order quantity of all periods for item 1 and item 3 of ant 1 will be copied into the solution generated by last ant. This novel concept represents the direct communication between ants in each cycle in order to further improve the solution constructed by individual ant separately.In this research work, ant based solution is constructed by using problem specific knowledge so that it does not violate constraint (3) and constraint (4), specified under mathematical modeling under Section2. Once the value of order quantity is generated for all components for all periods, it will be verified that whether the generated values of order quantity will meet the budget constraint and warehouse area constraint or not. If either of the constraints is not met, then the generated values of order quantity (solution) will be discarded. If both the budget constraint and the warehouse area constraint are met, then the solution is valid and it is stored for further processing.The ACO algorithm has been implemented on JAVA platform and the program can be run with Net Beans IDE.


### 3.1.2. Determination of ACO Parameter

To determine the optimum value for number of cycles, which is an important ACO algorithm parameter, ACO is run several times for one of the data sets, changing the number of cycles each time and recording the objective value function and CPU processing time. ACO run results are tabulated in Table1. From Figure 1 which shows the variation of objective function value with number of cycles in the ACO model, it is clear that objective function reaches optimum value at 600 cycles. After that increasing the number of cycles will only increase the CPU time without improving the objective function value.Table 1
Variation of objective function value and time of execution with number of cycles.
Number of cycles
Objective function value (Rs)
Time of execution (Sec)
50
8382155
5
100
8255916
9
200
8101560
20
300
8094680
32
500
8088430
41
600
8082500
54
700
8082435
62
800
8082440
75
1000
8082440
90Figure 1
Variation of objective function value with number of cycles: ACO model.


## 3.1.1. Ant Cycle Model

Ant cycle model is adopted in this research work. In this model, the trails are globally updated during each cycle by all ants. The amount of pheromone deposited by each ant is a function of the solution quality. As per the flowchart of proposed ant colony optimization algorithm, first the ant based initial solution construction method is executed. Next values of pheromone are set based on ant based initial solution construction method. Thus the values of pheromone are set accurately after the execution of ant based initial solution construction method. Followed by this, ant based solution construction method will be executed in each cycle which will utilize the values of pheromone in order to construct better solutions [11].This research work adopts ant cycle model where the trails are uniformly updated during each cycle by all ants. The solution quality is determined by the amount of pheromone deposited by each ant. The initial solution construction method is executed. Then the ant based solution construction method is executed which will employ the values of pheromone to construct better solutions. Each ant represents a solution to the optimization problem. In our case of application to inventory management, it represents the order quantity of different items in different periods so that total inventory cost is reduced.The novel concept is implemented for ant based solution construction in each cycle as follows. First each ant will construct the solution according to the probabilistic rule. After constructing the solution according to the probabilistic rule, ant 1 will discard its order quantity values of all periods for item 1 and item 3 and the values of order quantity of all periods for item 1 and item 3 of ant 2 will be copied into the solution generated by ant 1. Similarly, ant 2 will discard its order quantity values of all periods for item 1 and item 3 and the values of order quantity of all periods for item 1 and item 3 of ant 3 will be copied into the solution generated by ant 2. This discard-copy step is repeated for all ants. Last ant will discard its order quantity values of all periods for item 1 and item 3 and the values of order quantity of all periods for item 1 and item 3 of ant 1 will be copied into the solution generated by last ant. This novel concept represents the direct communication between ants in each cycle in order to further improve the solution constructed by individual ant separately.In this research work, ant based solution is constructed by using problem specific knowledge so that it does not violate constraint (3) and constraint (4), specified under mathematical modeling under Section2. Once the value of order quantity is generated for all components for all periods, it will be verified that whether the generated values of order quantity will meet the budget constraint and warehouse area constraint or not. If either of the constraints is not met, then the generated values of order quantity (solution) will be discarded. If both the budget constraint and the warehouse area constraint are met, then the solution is valid and it is stored for further processing.The ACO algorithm has been implemented on JAVA platform and the program can be run with Net Beans IDE.


## 3.1.2. Determination of ACO Parameter

To determine the optimum value for number of cycles, which is an important ACO algorithm parameter, ACO is run several times for one of the data sets, changing the number of cycles each time and recording the objective value function and CPU processing time. ACO run results are tabulated in Table1. From Figure 1 which shows the variation of objective function value with number of cycles in the ACO model, it is clear that objective function reaches optimum value at 600 cycles. After that increasing the number of cycles will only increase the CPU time without improving the objective function value.Table 1
Variation of objective function value and time of execution with number of cycles.
Number of cycles
Objective function value (Rs)
Time of execution (Sec)
50
8382155
5
100
8255916
9
200
8101560
20
300
8094680
32
500
8088430
41
600
8082500
54
700
8082435
62
800
8082440
75
1000
8082440
90Figure 1
Variation of objective function value with number of cycles: ACO model.


## 3.2. Genetic Algorithm

Genetic algorithm is a population based metaheuristic which can be utilized to solve complex optimization problems. John Holland invented genetic algorithm which mimics the principle of natural genetics in order to solve optimization problems. Genetic algorithm uses knowledge from previous generations to direct the future search. A number of chromosomes rather than a single chromosome can be used to carry out the search process. The optimum or near-optimum solution can be efficiently obtained using genetic algorithm [12].The GA implementation comprises the following important stages [13]:(i)
Setting GA parameters including the crossover probability (Pc), the mutation probability (Pm), population size (PS), and number of generations (NG).(ii)
Initializing the population with the size of PS randomly.(iii)
Evaluating the objective function.(iv)
Selecting individual for mating pool by roulette wheel selection method and using elitisms.(v)
Applying the crossover operation for each pair of chromosomes based on Pc.(vi)
Applying mutation operation for each chromosome based on Pm.(vii)
Replacing the current population by the resulting new population.(viii)
Continuing the cycle till stopping criterion is met, like the maximum number of cycles.The GA has been implemented on JAVA platform and the program can be run with Net Beans IDE. The following components of genetic algorithm have been implemented to solve the problem.


### 3.2.1. Chromosome Representation

The chromosomes represent the strings of the order quantities of the items in each period in this research work. For example, the chromosome representation with four items and three periods is shown in (14). In (14), the number of items and periods are represented by rows and columns, respectively.(14)CR=Q11Q12Q13Q21Q22Q23Q31Q32Q33Q41Q42Q43.


### 3.2.2. Initial Population Generation

Initial population of size PS is determined using Taguchi’s design which is explained in Section3.3.


### 3.2.3. Calculate Fitness of Each Chromosome

The process of calculating the fitness of each chromosome consists of the following two steps.(1)
Calculate the objective function of each chromosome.(2)
Calculate the fitness of each chromosome based on its objective function value using the following equation:(15)fx=11+Ox.In (15), fx is the fitness function of each chromosome and Ox is the objective function of each chromosome.


### 3.2.4. Selection Operator

The roulette wheel selection is used to select the chromosomes. In roulette wheel selection, chromosomes are selected into the mating pool according to their raw fitness. The chromosomes having high fitness have more chances to be selected. The probability ofith-selected chromosome is given by the following equation:(16)Pi=fi∑j=1nfj.In (16), f(i) and f(j) are the fitness of chromosomes i and j, respectively [12].


### 3.2.5. Crossover Operator

Single point crossover is used in this research work. The following equation (17) demonstrates the application of single point crossover on the two selected parents. Crossover rate is also optimized using Taguchi’s design in order to protect some of good chromosomes that are already present in the mating pool.(17)Parent1=Q11Q12Q13Q21Q22Q23Q31Q32Q33Q41Q42Q43,Offspring1=Q11Q12O13Q21Q22O23Q31Q32O33Q41Q42O43»Parent2=O11O12O13O21O22O23O31O32O33O41O42O43,Offspring2=O11O12Q13O21O22Q23O31O32Q33O41O42Q43.


### 3.2.6. Mutation Operator

A novel mutation operation, namely, “subtraction-addition” mutation, is used in this research work as follows.Step 1.
Any one period is selected randomly for each component.Step 2.
Mutation value is randomly generated.Step 3.
If selected period ist then mutation value is subtracted from period t and added to period T. If selected period t is the last period then T is the first period. If selected period t is not the last period then T=t+1 period.This novel mutation operation never violates constraint (3) of the total cost optimizing function which is stated in Section2.1.4. Again the mutation value is generated in a way that this novel mutation operation never violates constraint (4) too. Moreover, mutation rate is also optimized using Taguchi design.


### 3.2.7. Replacement Strategy

The chromosomes in old population are replaced by newly generated chromosomes which meet both the budget constraint (5) and warehouse area constraint (6) in each generation. The newly generated chromosomes which do not meet either constraint (5) or constraint (6) will be discarded.


### 3.2.8. Stopping Criteria

Genetic algorithm stops when it has completed maximum number of generations. Maximum number of generations can be increased or decreased based on the input size of the problem. This is another parameter of GA which is optimized using Taguchi experiment design.


## 3.2.1. Chromosome Representation

The chromosomes represent the strings of the order quantities of the items in each period in this research work. For example, the chromosome representation with four items and three periods is shown in (14). In (14), the number of items and periods are represented by rows and columns, respectively.(14)CR=Q11Q12Q13Q21Q22Q23Q31Q32Q33Q41Q42Q43.


## 3.2.2. Initial Population Generation

Initial population of size PS is determined using Taguchi’s design which is explained in Section3.3.


## 3.2.3. Calculate Fitness of Each Chromosome

The process of calculating the fitness of each chromosome consists of the following two steps.(1)
Calculate the objective function of each chromosome.(2)
Calculate the fitness of each chromosome based on its objective function value using the following equation:(15)fx=11+Ox.In (15), fx is the fitness function of each chromosome and Ox is the objective function of each chromosome.


## 3.2.4. Selection Operator

The roulette wheel selection is used to select the chromosomes. In roulette wheel selection, chromosomes are selected into the mating pool according to their raw fitness. The chromosomes having high fitness have more chances to be selected. The probability ofith-selected chromosome is given by the following equation:(16)Pi=fi∑j=1nfj.In (16), f(i) and f(j) are the fitness of chromosomes i and j, respectively [12].


## 3.2.5. Crossover Operator

Single point crossover is used in this research work. The following equation (17) demonstrates the application of single point crossover on the two selected parents. Crossover rate is also optimized using Taguchi’s design in order to protect some of good chromosomes that are already present in the mating pool.(17)Parent1=Q11Q12Q13Q21Q22Q23Q31Q32Q33Q41Q42Q43,Offspring1=Q11Q12O13Q21Q22O23Q31Q32O33Q41Q42O43»Parent2=O11O12O13O21O22O23O31O32O33O41O42O43,Offspring2=O11O12Q13O21O22Q23O31O32Q33O41O42Q43.


## 3.2.6. Mutation Operator

A novel mutation operation, namely, “subtraction-addition” mutation, is used in this research work as follows.Step 1.
Any one period is selected randomly for each component.Step 2.
Mutation value is randomly generated.Step 3.
If selected period ist then mutation value is subtracted from period t and added to period T. If selected period t is the last period then T is the first period. If selected period t is not the last period then T=t+1 period.This novel mutation operation never violates constraint (3) of the total cost optimizing function which is stated in Section2.1.4. Again the mutation value is generated in a way that this novel mutation operation never violates constraint (4) too. Moreover, mutation rate is also optimized using Taguchi design.


## 3.2.7. Replacement Strategy

The chromosomes in old population are replaced by newly generated chromosomes which meet both the budget constraint (5) and warehouse area constraint (6) in each generation. The newly generated chromosomes which do not meet either constraint (5) or constraint (6) will be discarded.


## 3.2.8. Stopping Criteria

Genetic algorithm stops when it has completed maximum number of generations. Maximum number of generations can be increased or decreased based on the input size of the problem. This is another parameter of GA which is optimized using Taguchi experiment design.


## 3.3. Parameter Calibration Using Taguchi Design of Experiments

Selection of optimum parameters apparently has significant effect on the effectiveness of the metaheuristic algorithms. When the number of factors affecting the performance of genetic algorithm is more, the full factorial experiment and exhaustive approach to investigate the effect of different parameters becomes increasingly complicated and impractical.In this work, the Taguchi method is used to tune the parameters of metaheuristic algorithm of multiobjective optimization, genetic algorithm. The parameters of a metaheuristic that are needed to be tuned act like controllable factors in the design of experiments (DOE). The aim is to find an optimal combination of the parameters such that the response (the fitness function) is optimized. The Taguchi method is a special case of the fractional factorial design in which some special orthogonal arrays are used [14]. Orthogonal arrays under Taguchi method would help to study a large number of decision variables with a limited number of experiments. Decision variables are divided into controllable and noise factors. Noise factors can not be controlled directly. It is also impractical and most of the time impossible to eliminate the noise factors [15]. Taguchi experimental design will help to reduce the effect of noise factors.Taguchi adopted the concept of signal to noise ratio to reduce the effect of noise factors in the experiment. The desired value or mean response value is represented by signal. The undesirable value or standard deviation is denoted by noise. The variation present in the response variable or the component of noise factor is represented byS/N ratio. Objective functions are classified into 3 types for design of experiment applications by Taguchi. They are “smaller the better,” the “larger the better,” and “the nominal is best.” Since almost all objective functions in inventory control systems are grouped in the “smaller the better” type, its corresponding S/N ratio is given by the following equation:(18)SN=-10log10⁡k2,where k is the objective function value.


### 3.3.1. Taguchi Method Implementation

Taguchi method implementation consists of five basic steps [14].(i)
Parameters with significant effects on the response are determined.(ii)
The parameter value is determined by the trial and error procedure so that good fitness value is obtained while implementing the experiments.(iii)
Taking into consideration the available degree of freedom (DOF), a suitable orthogonal array is selected which would specify the number of experiments to minimize the experimentation time and at the same time would be giving full weightage for all the combination of influencing factors.(iv)
Experiments are conducted based on obtained design.(v)
The results are recorded. TheS/N approach will give the basis for evaluation and analysis of the results.Setting GA parameters including the crossover probability (Pc), the mutation probability (Pm), population size (PS), and number of generations (NG) is very important in determining the efficiency of the metaheuristic algorithms like GA. Taguchi design which is based on fractional factorial experiments is used to set the GA parameters (see Table2). Minitab 15 is used to employ the Taguchi method. For 3 levels of 4 factors, L9 orthogonal array is suggested which will facilitate representative Fractional Factorial Experimentation. L9 orthogonal array lists the different combinations of factors at different levels at which the response value of experiments has to be determined. 3 items, 3 periods, and 3 price breaks data are selected for experimentation and response in terms of minimum total cost is tabulated as shown in Table 3. For different combinations of the factor levels, each example is solved three times and the mean response was used in the analysis. Figure 2 shows the main effect plot of S/N ratio for different parameter levels of the proposed algorithms. Based on “smaller the better” definition for the signal to noise ratio, optimum values of the GA parameters are as shown in the Table 4.Table 2
GA parameters and levels for Taguchi design.
Levels
GA parameters
Crossover probability Pc
Mutation probability Pm
Population size PS
Number of generations NG
1
0.7
0.15
4000
500
2
0.8
0.20
5000
600
3
0.9
0.30
6000
700Table 3
Experimental response for Taguchi design.
Crossover  probability Pc
Mutation probability Pm
Population size PS
Number of generations NG
Response value obj. function
S
/
N ratio
0.7
0.15
4000
500
8070423
-
138.138
0.7
0.20
5000
600
8071445
-
138.139
0.7
0.30
6000
700
8070334
-
138.138
0.8
0.15
5000
700
8070656
-
138.138
0.8
0.20
6000
500
8072345
-
138.140
0.8
0.30
4000
600
8074558
-
138.142
0.9
0.15
6000
600
8066732
-
138.134
0.9
0.20
4000
700
8065550
-
138.133
0.9
0.30
5000
500
8058597
-
138.125Table 4
Optimum GA parameters.
Optimum GA parameters
Crossover probability Pc
Mutation probability Pm
Population size PS
Number of generations NG
0.9
0.30
5000
500Figure 2
Main effect plot ofS/N ratios: Taguchi design.Minitab 15 is used to employ the Taguchi method. Under the menu options of Minitab, Stat-DOE-Taguchi Design-Create Design is selected. For 3 levels of 4 factors, L9 orthogonal array is suggested which will facilitate representative Fractional Factorial Experimentation. L9 orthogonal array lists the different combinations of factors at different levels at which the response values of experiments have to be determined. 3 items, 3 periods, and 3 price breaks data are selected for experimentation.For different combinations of the factor levels, each example is solved three times and the mean response was used in the analysis. Figure2 shows the main effect plot of S/N ratio and Figure 3 shows the main effect plot of means for different parameter levels of the proposed algorithms.Figure 3
Main effect plot of means: Taguchi design.TheS/N ratio indicates the amount of variation present in the response variable, and the aim is to maximize it. Pc value of 0.9, Pm value of 0.3, PS value of 5000, and NG value of 500 yield the maximum value of S/N ratio as can be seen from Figure 2. This can also be verified from Figure 3 where the above values give the best or the lowest value of objective function. Optimal parameter values of the algorithms are shown in Table 4.


## 3.3.1. Taguchi Method Implementation

Taguchi method implementation consists of five basic steps [14].(i)
Parameters with significant effects on the response are determined.(ii)
The parameter value is determined by the trial and error procedure so that good fitness value is obtained while implementing the experiments.(iii)
Taking into consideration the available degree of freedom (DOF), a suitable orthogonal array is selected which would specify the number of experiments to minimize the experimentation time and at the same time would be giving full weightage for all the combination of influencing factors.(iv)
Experiments are conducted based on obtained design.(v)
The results are recorded. TheS/N approach will give the basis for evaluation and analysis of the results.Setting GA parameters including the crossover probability (Pc), the mutation probability (Pm), population size (PS), and number of generations (NG) is very important in determining the efficiency of the metaheuristic algorithms like GA. Taguchi design which is based on fractional factorial experiments is used to set the GA parameters (see Table2). Minitab 15 is used to employ the Taguchi method. For 3 levels of 4 factors, L9 orthogonal array is suggested which will facilitate representative Fractional Factorial Experimentation. L9 orthogonal array lists the different combinations of factors at different levels at which the response value of experiments has to be determined. 3 items, 3 periods, and 3 price breaks data are selected for experimentation and response in terms of minimum total cost is tabulated as shown in Table 3. For different combinations of the factor levels, each example is solved three times and the mean response was used in the analysis. Figure 2 shows the main effect plot of S/N ratio for different parameter levels of the proposed algorithms. Based on “smaller the better” definition for the signal to noise ratio, optimum values of the GA parameters are as shown in the Table 4.Table 2
GA parameters and levels for Taguchi design.
Levels
GA parameters
Crossover probability Pc
Mutation probability Pm
Population size PS
Number of generations NG
1
0.7
0.15
4000
500
2
0.8
0.20
5000
600
3
0.9
0.30
6000
700Table 3
Experimental response for Taguchi design.
Crossover  probability Pc
Mutation probability Pm
Population size PS
Number of generations NG
Response value obj. function
S
/
N ratio
0.7
0.15
4000
500
8070423
-
138.138
0.7
0.20
5000
600
8071445
-
138.139
0.7
0.30
6000
700
8070334
-
138.138
0.8
0.15
5000
700
8070656
-
138.138
0.8
0.20
6000
500
8072345
-
138.140
0.8
0.30
4000
600
8074558
-
138.142
0.9
0.15
6000
600
8066732
-
138.134
0.9
0.20
4000
700
8065550
-
138.133
0.9
0.30
5000
500
8058597
-
138.125Table 4
Optimum GA parameters.
Optimum GA parameters
Crossover probability Pc
Mutation probability Pm
Population size PS
Number of generations NG
0.9
0.30
5000
500Figure 2
Main effect plot ofS/N ratios: Taguchi design.Minitab 15 is used to employ the Taguchi method. Under the menu options of Minitab, Stat-DOE-Taguchi Design-Create Design is selected. For 3 levels of 4 factors, L9 orthogonal array is suggested which will facilitate representative Fractional Factorial Experimentation. L9 orthogonal array lists the different combinations of factors at different levels at which the response values of experiments have to be determined. 3 items, 3 periods, and 3 price breaks data are selected for experimentation.For different combinations of the factor levels, each example is solved three times and the mean response was used in the analysis. Figure2 shows the main effect plot of S/N ratio and Figure 3 shows the main effect plot of means for different parameter levels of the proposed algorithms.Figure 3
Main effect plot of means: Taguchi design.TheS/N ratio indicates the amount of variation present in the response variable, and the aim is to maximize it. Pc value of 0.9, Pm value of 0.3, PS value of 5000, and NG value of 500 yield the maximum value of S/N ratio as can be seen from Figure 2. This can also be verified from Figure 3 where the above values give the best or the lowest value of objective function. Optimal parameter values of the algorithms are shown in Table 4.


## 4. Performance Evaluation and Comparison

Figure4 highlights performance measures that are adopted in the research work for evaluating effectiveness of GA and ACO models and for their comparison. In order to evaluate the performance of the developed models in suggesting the most economical periodical lot sizing, taking into consideration dynamic demand and ordering, holding, and purchase cost, 8 sample data sets are developed. The data input pertains to the real industry scenario of valve manufacturing company. Data for six most important valve configurations is studied in depth and results obtained by applying both algorithms are tabulated. Performance of two algorithms is evaluated based on 4 important parameters of minimum total cost or objective function value, CPU execution time, spread between best and worst solution in the program run, and sensitivity analysis of program parameters.Figure 4
ACO and GA performance measures.Each data set can be conveniently represented as data set number (s-t-p), where s is the number of items for procurement, t is the period, and p is the number of price breaks for discount. Typical input data for the lot sizing problem is shown in Table 5 which is for data set number 8 (s=6,t=3,p=3), where six-item procurement is considered for 3 time horizons with 3 price breaks.Table 5
Input data for lot sizing problem: data set number 8.
Data set number 8 (6-3-3)
s
=
6, t=3, p=3
s
6
D61
425
P11
6025
R1
20
t
3
D62
440
P12
5422.5
R2
40
p
3
D63
460
P13
5121.25
R3
85
D11
40
K1
25000
P21
5100
R4
30
D12
70
K2
20000
P22
4590
R5
190
D13
50
K3
9000
P23
4335
R6
130
D21
130
K4
5000
P32
4346.25
B1
10000000
D22
110
K5
7000
P33
4117.5
B2
10000000
D23
140
K6
7000
P41
3550
B3
10000000
D31
300
H1
600
P42
3195
m1
0.25
D32
260
H2
300
P43
3017
m2
0.2
D33
280
H3
425
P51
4500
m3
0.15
D41
80
H4
200
P52
4275
m4
0.1
D42
92
H5
150
P53
4050
m5
0.1
D43
98
H6
300
P61
7000
m6
0.15
D51
600
P62
6650
M
600
D52
665
P63
6300
D53
620ACO and GA were run 4 times for each problem configuration. For each run, the order quantities of every item for each period are tabulated. Minimum total cost or the value of objective function which is an important performance measure is listed. The CPU execution time is also noted for each run. The spread between best and worst solution is also marked for each run. The lesser the spread, the better the algorithm. Table6 shows the results for ACO run and Table 7 shows the results for GA run.Table 6
ACO performance results.
Problem type
Run number
Item number
Order qty. (number) period 1
Orderqty. (number) period 2
Order qty. (number) period 3
Best sol. total cost (Rs)
Spread bet. best and worst solution
CPU time (Sec)
s
=
1 t=3 p=2
1
1
51
59
50
1159882
144660
12
2
1
51
59
50
1159882
144660
12
3
1
51
59
50
1159882
144660
12
4
1
51
59
50
1159882
144660
12
s
=
1 t=3 p=3
1
1
40
120
0
1135500
121445
12
2
1
40
120
0
1135500
121445
12
3
1
40
120
0
1135500
121445
12
4
1
40
120
0
1135500
121445
12
s
=
2 t=3 p=2
1
1
51
59
50
3215982
336700
20
2
130
110
140
2
1
51
59
50
3215982
336700
20
2
130
110
140
3
1
51
59
50
3215982
336700
20
2
130
110
140
4
1
51
59
50
3215982
336700
20
2
130
110
140
s
=
2 t=3 p=3
1
1
40
120
0
3494700
316015
22
2
130
110
140
2
1
40
120
0
3494700
316015
22
2
130
110
140
3
1
40
120
0
3494700
316015
22
2
130
110
140
4
1
40
120
0
3494700
316015
22
2
130
110
140
s
=
3 t=3 p=3
1
1
42
118
131
7015841
493352
40
2
131
120
129
3
361
208
271
2
1
40
120
0
7010875
497825
40
2
134
108
138
3
372
468
00
3
1
43
117
0
7017807.
485675
40
2
138
102
140
3
372
468
0
4
1
41
119
0
7021825
491372
40
2
130
110
140
3
377
426
377
S
=
4 t=3 p=3
1
1
41
103
16
8070056.25
487646.25
62
2
141
109
130
3
375
214
251
4
270
0
0
2
1
60
53
47
8072996.25
469560.0
62
2
139
103
138
3
361
414
65
4
270
0
0
3
1
51
109
0
8082756.25
485410.0
62
2
144
112
124
3
377
401
62
4
270
0
0
4
1
54
106
0
8073832.5
489540.0
62
2
131
116
133
3
380
402
58
4
270
0
0
s
=
5 t=3 p=3
1
1
160
0
0
15718292
755940.25
75
2
380
0
0
3
361
244
235
4
270
0
0
5
756
654
475
2
1
51
109
0
16447257
720387.5
75
2
159
211
10
3
370
233
237
4
270
0
0
5
671
756
458
3
1
157
3
0
15726613
745871.25
75
2
380
0
0
3
365
211
264
4
110
144
16
5
815
504
566
4
1
53
107
0
15702395
766744
75
2
153
227
0
3
307
253
280
4
270
0
0
5
663
711
511
s
=
6 t=3 p=3
1
1
46
114
0
25321202.5
645500.0
80
2
153
210
17
3
377
411
52
4
260
6
4
5
628
721
536
6
480
413
432
2
1
72
55
33
25331901
622969
80
2
160
220
11
3
377
380
83
4
268
2
0
5
635
648
602
6
444
425
456
3
1
52
108
0
25324673
651288.75
80
2
203
161
16
3
364
262
214
4
104
146
20
5
616
680
589
6
429
440
456
4
1
42
118
0
25301616
680659
80
2
201
179
0
3
366
244
230
4
270
0
0
5
637
658
590
6
426
456
443Table 7
GA performance results.
Problem type
Run number
Item number
Order qty. (number) period 1
Order qty. (number) period 2
Order qty. (number) period 3
Best sol. total cost (Rs)
Spread bet. best and worst solution
CPU time (Sec)
s
=
1 t=3 p=2
1
1
51
59
50
1159882
144660
11
2
1
51
59
50
1159882
144660
11
3
1
51
59
50
1159882
144660
11
4
1
51
59
50
1159882
144660
11
s
=
1 t=3 p=3
1
1
40
120
0
1135500
121445
11
2
1
40
120
0
1135500
121445
11
3
1
40
120
0
1135500
121445
11
4
1
40
120
0
1135500
121445
11
s
=
2 t=3 p=2
1
1
51
59
50
3215982
336700
12
2
130
110
140
2
1
51
59
50
3215982
336700
12
2
130
110
140
3
1
51
59
50
3215982
336700
12
2
130
110
140
4
1
51
59
50
3215982
336700
12
2
130
110
140
s
=
2 t=3 p=3
1
1
40
120
0
3494700
316015
12
2
130
110
140
2
1
40
120
0
3494700
316015
12
2
130
110
140
3
1
40
120
0
3494700
316015
12
2
130
110
140
4
1
40
120
0
3494700
316015
12
2
130
110
140
s
=
3 t=3 p=3
1
1
42
118
131
7015841
493352
12
2
131
120
129
3
361
208
271
2
1
40
120
0
7010875
497825
12
2
134
108
138
3
372
468
00
3
1
43
117
0
7017807.
485675
12
2
138
102
140
3
372
468
0
4
1
41
119
0
7021825
491372
12
2
130
110
140
3
377
426
377
S
=
4 t=3 p=3
1
1
51
109
0
8090702
510110
12
2
131
109
140
3
367
227
246
4
80
106
84
2
1
41
119
0
8058597
514989
12
2
140
106
164
3
368
472
0
4
81
96
93
3
1
40
120
0
8067023
514867
12
2
130
113
137
3
367
228
245
4
80
119
79
4
1
40
120
0
8056263
498388
12
2
120
112
129
3
363
221
256
4
83
89
98
s
=
5 t=3 p=3
1
1
57
103
0
15764097
802316
13
2
168
212
0
3
407
395
38
4
270
0
0
5
634
742
509
2
1
45
115
0
15659722
828462
13
2
153
227
0
3
367
208
265
4
120
150
0
5
674
623
588
3
1
47
113
0
15736621
785438
13
2
164
216
0
3
312
254
274
4
270
0
0
5
828
446
611
4
1
48
112
0
15648478
790905
13
2
165
215
0
3
369
407
64
4
270
0
0
5
600
665
620
s
=
6 t=3 p=3
1
1
105
55
0
25260387
675846
13
2
170
210
0
3
364
209
267
4
270
0
0
5
620
651
614
6
426
442
457
2
1
54
106
0
2
159
221
0
25287605
692205
18
3
370
214
256
4
149
93
28
5
603
662
620
6
425
440
460
3
1
57
103
0
25307731.25
699552.5
13
2
157
223
0
3
430
387
23
4
150
120
0
5
617
650
618
6
428
438
459
4
1
81
79
0
25305686
658103
13
2
167
213
0
3
373
223
244
4
118
137
15
5
602
665
618
6
428
437
460For data sets 1 and 2 with only one component, where the number of variables is minimum, we see that all the four runs of GA and ACO return the same results for order quantities and objective function value. Our ACO and GA prescribed the same optimal solution each time and the same solution that was obtained via exact methods using the LINDO 14 Linear Integer Programming software. Since this is a case of simple problem which involves small number of variables, the solution obtained using ACO and GA is the exact or the most optimal one.Solving the above problem using other optimization software like GAMS or LINDO also gives the same result. This fact validates model. Having verified its ability to optimize the lot sizing, ACO and GA can now be used to solve the complex problems having more number of decision variables.However, exact methods can not be used for large size problem because of very large number of variables involved, the program becomes unstable or takes very long time for execution.For lot sizing decisions for 3 components and above, the number of decision variables increases and we can find different solutions in different ACO and GA runs. This is evident from Tables6 and 7. Four problems are considered with the same number of periods and price breaks but different number of items.One-way ANOVA analysis is carried out to compare the performance of GA and ACO model.Table8 lists the result of ANOVA analysis to compare the solution methodologies statistically on the basis of best solution cost or best fitness value and CPU execution time. The output indicates that, at confidence level 95%, the two algorithms have no statistically significant differences in the best fitness value between their means as P value is 0.987>0.05. But when CPU execution time is concerned, two algorithms are significantly different as P value is less than 0.05. This is based on the criteria that the null hypothesis that the population means are equal is accepted if P value is greater than 0.05, otherwise rejected.Table 8
ANOVA analysis to compare ACO & GA.
Source
DF
SS
MS
F
P value
Parameter
BFV
CPT
BFV
CPT
BFV
CPT
BFV
CPT
BFV
CPT
Solution methodology
1
1
18643854306
3306
18643854306
3306
0
8.26
0.987
0.012
Error
14
14
1.0069
7
E
+
15
5604
7.1926
7
E
+
13
400
Total
15
15
BFV: best fitness value. CPT: CPU processing time.From Tables6 and 7, it can be observed that GA is exceptionally better than ACO based on the performance parameter of CPU execution time. When the problem is simple with smaller number of variables, there is no much variance in the execution time. But as the problem becomes more complex and the number of variables increases, there is a significant difference and a trend of lesser CPU time in the case of GA model. This demonstrates the higher efficiency of GA compared to ACO algorithm in terms of CPU execution time for solving our mathematical model of inventory cost optimization. This can be attributed to the simplicity of GA compared to ACO in the algorithm formulation which is reflected in the lesser CPU time.


### 4.1. Sensitivity Analysis

Sensitivity analysis is the key to performance appraisal of the mathematical models and their solutions [16]. Sensitivity analysis establishes and verifies the relationship between problem parameters and thereby validating the model. In this section, the effect of variation in problem parameters is investigated on multi-period lot sizing decision from our algorithm. The influence of variation in ordering cost, holding cost, and purchase cost in terms of discounts and price breaks on the order quantity and minimum cost is studied for the proposed algorithm.Simple problem with data set number 1 (s = 1, t=3, p=2), where single item procurement is considered over 3 time horizons and 2 price breaks, is used to study the effect of various parameters on procurement lot size and total cost. It was noticed that both ACO and GA run gave the same results as the number of variables was limited.Table9 and Figure 5 clearly show that as the ordering cost is increased, the optimum number of orders decreases. The algorithm suggested order quantities for 3 periods when the ordering cost was less than around Rs 35940/-. After that, the lot sizing suggested to buy the quantities only at 2 periods reducing the number of annual orders to 2. When the ordering cost is still higher at Rs 111180/-, the algorithm suggests only one lot to be bought at the year beginning on the basis of ordering cost-holding cost tradeoff. The saving in ordering cost was more than the increase in holding cost expenses.Table 9
Ordering cost-sensitivity analysis.
Run number
Ordering cost (Rs)
Order qty. (number) period 1
Order qty. (number) period 2
Order qty. (number) Period 3
Best sol. total cost (Rs)
1
25000
51
59
50
1159882
2
35000
51
59
50
1189882
3
35940
51
109
0
1192700
4
111180
160
0
0
1343180Figure 5
Ordering cost: sensitivity analysis.Table10 and Figure 6 show the different periodic lot size by running the GA and ACO by varying the holding cost. For the present problem, at holding cost of Rs 1020/-, the order quantities are 51, 59, and 50 for the three periods. If the holding cost is increased, a point will be reached at around Rs 1100/-, where the GA and ACO algorithm will suggest the JIT (Just in Time) procurement of quantities 40, 70, and 50, respectively, for the 3 periods because, at this cost, there can not be any tradeoff between holding and ordering cost.Table 10
Holding cost -sensitivity analysis.
Run number
Ordering cost (Rs)
Order qty. (number) period 1
Order qty. (number) period 2
Order qty. (number) period 3
Best sol. total cost (Rs)
1
1020
51
59
50
1159882
2
1100
40
70
50
1171912
3
800
51
109
0
1126600
4
225
160
0
0
1010550Figure 6
Holding cost-sensitivity analysis.It will be interesting to note that as the holding cost is reduced to around Rs 800/-, the first tradeoff point between ordering and holding cost is reached. So it becomes more advantageous to reduce one order in spite of higher carrying expenses due to larger inventory. So optimized order quantities are obtained as 51, 109, and 0.When the holding cost is further reduced, at around Rs 220, it becomes less expensive to buy all the annual demand of 160 at the beginning of the year, due to saving in ordering cost.Table11 records the sensitivity of the order quantity decisions on the discounts and price breaks.Table 11
Sensitivity to discounts and price breaks.
Number
Discount/price break
Discount/price break 2
Order qty. (number) period 1
Order qty. (number) period 2
Order qty. (number) period 3
Best sol. total cost (Rs)
1
5% at & above 50
NIL
51
59
50
1159882
2
5% at & above 60
10% at & above 100
40
120
0
1135500
3
25% at & above 150
Nil
160
0
0
1064200In scenario number 1, all the order quantities are so selected that they are above 50 to see that they are eligible for the 5% discount at the first price break. In scenario 3, there is single discount of 25% above price break quantity of 150. To avail this huge discount, the order quantity of 160 is suggested in the first period itself to ensure that the order quantity stays at above 150.It is evident from above sensitivity analysis that procurement lot sizing model proposed is sensitive to the variation in problem parameters. The computational results suggest that the proposed model captures all realistic constraints in multi-period procurement lot sizing decision-making process and analyzes tradeoffs in cost objectives. Optimal procurement lot size is obtained by striking best tradeoffs among multiple cost objectives.Smaller lot size reduces inventory holding cost but increases purchasing cost and ordering cost due to lack of economy of scale. Larger lot size reduces purchasing cost and leads to higher inventory holding cost. All these are reflected in our sensitivity analysis.


## 4.1. Sensitivity Analysis

Sensitivity analysis is the key to performance appraisal of the mathematical models and their solutions [16]. Sensitivity analysis establishes and verifies the relationship between problem parameters and thereby validating the model. In this section, the effect of variation in problem parameters is investigated on multi-period lot sizing decision from our algorithm. The influence of variation in ordering cost, holding cost, and purchase cost in terms of discounts and price breaks on the order quantity and minimum cost is studied for the proposed algorithm.Simple problem with data set number 1 (s = 1, t=3, p=2), where single item procurement is considered over 3 time horizons and 2 price breaks, is used to study the effect of various parameters on procurement lot size and total cost. It was noticed that both ACO and GA run gave the same results as the number of variables was limited.Table9 and Figure 5 clearly show that as the ordering cost is increased, the optimum number of orders decreases. The algorithm suggested order quantities for 3 periods when the ordering cost was less than around Rs 35940/-. After that, the lot sizing suggested to buy the quantities only at 2 periods reducing the number of annual orders to 2. When the ordering cost is still higher at Rs 111180/-, the algorithm suggests only one lot to be bought at the year beginning on the basis of ordering cost-holding cost tradeoff. The saving in ordering cost was more than the increase in holding cost expenses.Table 9
Ordering cost-sensitivity analysis.
Run number
Ordering cost (Rs)
Order qty. (number) period 1
Order qty. (number) period 2
Order qty. (number) Period 3
Best sol. total cost (Rs)
1
25000
51
59
50
1159882
2
35000
51
59
50
1189882
3
35940
51
109
0
1192700
4
111180
160
0
0
1343180Figure 5
Ordering cost: sensitivity analysis.Table10 and Figure 6 show the different periodic lot size by running the GA and ACO by varying the holding cost. For the present problem, at holding cost of Rs 1020/-, the order quantities are 51, 59, and 50 for the three periods. If the holding cost is increased, a point will be reached at around Rs 1100/-, where the GA and ACO algorithm will suggest the JIT (Just in Time) procurement of quantities 40, 70, and 50, respectively, for the 3 periods because, at this cost, there can not be any tradeoff between holding and ordering cost.Table 10
Holding cost -sensitivity analysis.
Run number
Ordering cost (Rs)
Order qty. (number) period 1
Order qty. (number) period 2
Order qty. (number) period 3
Best sol. total cost (Rs)
1
1020
51
59
50
1159882
2
1100
40
70
50
1171912
3
800
51
109
0
1126600
4
225
160
0
0
1010550Figure 6
Holding cost-sensitivity analysis.It will be interesting to note that as the holding cost is reduced to around Rs 800/-, the first tradeoff point between ordering and holding cost is reached. So it becomes more advantageous to reduce one order in spite of higher carrying expenses due to larger inventory. So optimized order quantities are obtained as 51, 109, and 0.When the holding cost is further reduced, at around Rs 220, it becomes less expensive to buy all the annual demand of 160 at the beginning of the year, due to saving in ordering cost.Table11 records the sensitivity of the order quantity decisions on the discounts and price breaks.Table 11
Sensitivity to discounts and price breaks.
Number
Discount/price break
Discount/price break 2
Order qty. (number) period 1
Order qty. (number) period 2
Order qty. (number) period 3
Best sol. total cost (Rs)
1
5% at & above 50
NIL
51
59
50
1159882
2
5% at & above 60
10% at & above 100
40
120
0
1135500
3
25% at & above 150
Nil
160
0
0
1064200In scenario number 1, all the order quantities are so selected that they are above 50 to see that they are eligible for the 5% discount at the first price break. In scenario 3, there is single discount of 25% above price break quantity of 150. To avail this huge discount, the order quantity of 160 is suggested in the first period itself to ensure that the order quantity stays at above 150.It is evident from above sensitivity analysis that procurement lot sizing model proposed is sensitive to the variation in problem parameters. The computational results suggest that the proposed model captures all realistic constraints in multi-period procurement lot sizing decision-making process and analyzes tradeoffs in cost objectives. Optimal procurement lot size is obtained by striking best tradeoffs among multiple cost objectives.Smaller lot size reduces inventory holding cost but increases purchasing cost and ordering cost due to lack of economy of scale. Larger lot size reduces purchasing cost and leads to higher inventory holding cost. All these are reflected in our sensitivity analysis.


## 5. Conclusions and Directions of Future Research

In this research work, a multi-item multi-period procurement lot sizing inventory control problem was investigated to obtain the optimum order quantity so that total inventory cost is minimized. A mixed binary integer mathematical programming model was developed to represent the problem domain. More practical considerations with budget and storage constraints, price discounts, order size limitations, and variable periodic demand rate were incorporated. AI tools were adopted to solve these multiobjective optimization scenarios of inventory management. Two population based metaheuristics, variation of ACO and advanced GA, were utilized to solve the proposed NP-hard problem. For simpler problems with less number of variables, the solutions obtained were the same as solutions from exact algorithms run by GAMS software. This fact helps to validate the proposed model.Taguchi method was implemented to calibrate the parameters of the metaheuristics. Using data from real time industrial scenario, the performances of the two metaheuristic methodologies were compared based on two measures: the objective function value and the required CPU time to obtain the solution. There is no difference in the performance of both ACO and GA methods in terms of best objective function value achievable for inventory planning, but the GA method is far more fast with lesser CPU processing time.Future research direction may be towards exploring other metaheuristic search algorithms such as simulated annealing, which may be employed for the optimization of inventory planning and procurement lot sizing. Uncertainty in the estimation of the different variables like carrying cost, ordering cost, and others can be modeled by the fuzziness to take care of their stochastic nature which will give a different approach to problem solution.



---

*Source: 3601217-2017-12-20.xml*