A Particle Swarm Optimization hyper-heuristic for the Dynamic Vehicle Routing Problem
aa r X i v : . [ c s . N E ] J un A Particle Swarm Optimization hyper-heuristicfor the Dynamic Vehicle Routing Problem
Micha l Okulewicz
Faculty of Mathematics and Information ScienceWarsaw University of Technology, Warsaw, Poland
Jacek Ma´ndziuk
Faculty of Mathematics and Information ScienceWarsaw University of Technology, Warsaw, Poland
Abstract
This paper presents a method for choosing a Particle Swarm Opti-mization based optimizer for the Dynamic Vehicle Routing Problemon the basis of the initially available data of a given problem instance.The optimization algorithm is chosen on the basis of a prediction madeby a linear model trained on that data and the relative results obtainedby the optimization algorithms. The achieved results suggest that sucha model can be used in a hyper-heuristic approach as it improved theaverage results, obtained on the set of benchmark instances, by choos-ing the appropriate algorithm in 82% of significant cases. Two leadingmulti-swarm Particle Swarm Optimization based algorithms for solvingthe Dynamic Vehicle Routing Problem are used as the basic optimiza-tion algorithms: Khouadjia’s et al. Multi–Environmental Multi–SwarmOptimizer and authors’ 2–Phase Multiswarm Particle Swarm Optimiza-tion.
Keywords:
Dynamic Vehicle Routing Problem,Particle Swarm Optimization,Hyper-heuristic
1. Introduction
Dynamic transportation problems have been considered in the liter-ature by Psaraftis [21, 22] since 1980. After the introduction of a setbenchmarks by Kilby [15] and Montemanni [16], several meta-heuristicbased algorithms have been developed in order to solve the problem, in- cluding Ant Colony Optimization (ACS) [7,16], Genetic Algorithm (GA)[8, 10], Tabu Search (TS) [10], and Particle Swarm Optimization (PSO)[14, 18].Although some of the works [8,10,13,17] mention the features of a spa-tial distribution of the requests of a given DVRP instance (describing itas spatially uniform or clustered), none of those methods directly use theinformation about the known requests location and volume. An initialnon-parametric approach for using the information about the availablerequests volumes in order to generate artificial requests (to account forthe unknown ones) has been presented by the authors [19].This paper proposes a hyper-heuristic method based on Multi–Envi-ronmental Multi–Swarm Optimizer (MEMSO) [12–14] and 2–Phase Mul-tiswarm Particle Swarm Optimization (2MPSO) [17–19] algorithms. Thehyper-heuristic uses the statistical data about the initially known set ofrequests in the given Dynamic Vehicle Routing Problem (DVRP).The rest of the paper is organized as follows. Section 2 defines the DVRPsolved in this paper. Section 3 introduces PSO and algorithms MEMSOand 2MPSO (both based on the PSO) used for optimizing DVRP. Section4 presents the hyper-heuristic approach for solving the DVRP. Section5 gives experimental setup and results obtained by the method. Finally,Section 6 concludes the paper.
2. Dynamic Vehicle Routing Problem
DVRP is a dynamic version of the generalization of a Traveling Sales-man Problem (TSP), called Vehicle Routing Problem (VRP). In theVRP the goal is to optimize a total route for a fleet of vehicles witha limited capacity. VRP has been introduced as a
Truck dispatchingproblem in 1959 by Dantzig and Ramser [6]. After the technologicaladvancement of the vehicle tracking devices and development of the Ge-ographical Information Systems a notion of a DVRP has been reintro-duced by Psaraftis [22]. The problem has attracted more attention aftera set of static VRP benchmarks of Christofides [2], Fisher [9] and Tail-lard [24] has been customized for the DVRP by Kilby et al. in 1998 [15]and Montemanni et al. in 2005 [16].In this paper a most common variant of the DVRP is solved [20],sometimes referred to as a VRP with Dynamic Requests (VRPwDR)[13]. In this variant a homogeneous fleet of vehicles (identical capacity c ∈ R and speed sp ∈ R ) is considered. There is also an additional con-straint, that the vehicle may operate only during a working day definedby the opening hours of its depot. During that working day a fleet of m vehicles must serve (visit) a set of n requests. Each request is defined by Particle Swarm Optimization hyper-heuristic for the Dynamic Vehicle Routing Problem a location l i ∈ R , an amount of cargo s i (0 ≤ s i ≤ c ) to be delivered andan amount of time u i ∈ R it takes to provide a service at that location(unload that cargo). The dynamic nature of that optimization problemcomes from the fact that new request may arrive during the working day.The period of time during which new requests are accepted is limited bya cut–off time parameter T co . In all the benchmarks considered in thispaper there is one depot for all the vehicles and T co is equal to a half ofthe working day.In the mathematical sense, the DVRP is a problem of finding a set oflocation permutations resulting in a shortest total path length under thetime and capacity constrains for each of the permutations in which all thelocations are visited exactly once. The exact mathematical formulationof the problem can be found in [8, 18, 19].
3. MEMSO and 2MPSO Algorithms
In this section two most successful approaches to solving DVRP:Khouadjia’s et al. MEMSO and authors’ 2MPSO are presented. Thosetwo methods are the base algorithms for the hyper–heuristic approachproposed in this paper.Both, the MEMSO and the 2MPSO, use PSO as their base meta-heuristic and 2-OPT [5] as a route optimization heuristic. In bothmethods the working day is divided into discrete number of time slices,with the instance of the DVRP problem ”frozen” within each time slice.Therefore, each method solves a series of dependent static VRP instancesduring the optimization process. The difference between PSO applica-tions, encoding of the problem and knowledge transfer in the MEMSOand the 2MPSO methods, together with a brief description of the PSOand 2-OPT algorithms, are presented and discussed in this section.
PSO algorithm is an iterative population based continuous optimiza-tion meta-heuristic approach utilizing the concept of Swarm Intelligence.The algorithm has been introduced by Kennedy and Eberhart in 1995[11] and has been further developed and studied by other researchers[3, 4, 23].During the optimization process PSO maintains a set of fitness func-tion solutions (called particles). Each particle has its own location x ∈ R n (an n − dimensional fitness function solution proposal), velocity v (a solution change vector), a set of neighbors N (particles which solu- tions it can observe), a memory of the best observed solution x ( BEST ) N and a memory of the best visited solution x ( BEST ) i .In each iteration t the location vector x and velocity vector v of i thparticle are changed in the following way: x i,t = x i,t − + v i,t − (1) v i,t = ωv i,t − + u c ( x ( BEST ) i − x i,t ) + u c ( x ( BEST ) N − x i,t ) (2)Where ω denotes an inertia factor, c and c are personal and globalattraction factors, u and u follow the uniform n − dimensional distri-bution on [0 , n . Figure 1.
Depiction of a 2–OPT algorithm optimization process over a sampledirected cycle. While considering edges BE and CF the algorithm observed thatchanging them to BC and EF results in a shorter route. After swapping the endsof the considered edges a final step of reversing the direction of the E-D-C path isperformed.
Particle Swarm Optimization hyper-heuristic for the Dynamic Vehicle Routing Problem MEMSO [14] algorithm uses PSO to optimize division of the requestsamong the vehicles. The vehicles’ routes within those divided sets arecreated by a greedy insertion and optimized by the 2-OPT algorithm.The fitness function value is the total length of those routes.MEMSO uses a discrete encoding of the requests division. The so-lution is an integer vector representing the requests and the values inthe vector are the vehicles’ identifiers. Therefore, the PSO algorithm ischanged in such a way that a velocity is a vector in { , , . . . , m } n andaddition in eq. (1) is performed in Z mn space (where m is a number ofvehicles and n a number of requests).The knowledge about the current solution is transferred between sub-sequent time slices by adapting the whole population. For each of theparticles the already served and decisively assigned requests are blockedfrom being changed and new requests are inserted in a greedy way intothe solution vectors. k beingthe parameter of the method. Therefore, the PSO particle for the firstphase optimizer is a sequence of requests clusters centers flattened intoa vector in R k ˆ m space (where ˆ m is the estimated number of vehiclesnecessary to serve the requests). Particle in the PSO instance optimizingthe route of the vehicle is a sequence of requests ranks. Therefore, itis a vector in R n i (where n i is the number of requests assigned to i thvehicle).The knowledge is transferred between subsequent time slices in a formof cluster centers vector generating a division of requests and a set of requests rank vectors for the routes imposed by requests division. Thetransferred solution is expanded by a new random cluster center if morevehicles seem to be necessary and the initial ranks of the new requestsare also initialized at random.
4. Hyper-heuristic Approach
Figure 2.
The activity diagram presenting a single run of a 2MPSO and MEMSObased hyper-heuristic.
Although 2MPSO outperforms MEMSO by 2.22% on average on theKilby’s [15] and Montemanni’s [16] sets of benchmarks (see Table 2),it might be limited by its clustering approach (or needs an excessivelylarge k ) for some particular DVRP tasks. Particle Swarm Optimization hyper-heuristic for the Dynamic Vehicle Routing Problem For that reason, the authors propose a hyper-heuristic approach [1] inwhich a statistical model might be incrementally trained for choosing thealgorithm, which seems to be most suitable for a given DVRP instance.A single run of such hyper-heuristic (solving a single DVRP task) isdepicted in Fig. 2. Please observe in the activity diagram, that only thechosen algorithm is used to provide an actual output to some externaldecision support system. The other algorithm is run only to gather thedata and its result is used to tune the prediction model. Please alsonote, that the choice of the optimization algorithm is done only at thebeginning of the optimization process. The prediction model is trainedon the results gathered through a subsequent runs on different DVRPtasks.The economic cost of running such a hyper-heuristic, in comparisonwith a single algorithm optimization, would be slightly more than dou-bled. The doubling comes from the fact, that it is necessary to makea similar amount of computations by two algorithms in order to get com-parable results for the performance prediction model. The additionaloverhead is a result of computations needed for getting statistics fromthe set of requests and creating a prediction model over the already com-puted cases. Although the economic cost of proposed approach woulddefinitely be larger than that of a single algorithm, the overhead for thecomputations necessary for getting solutions during daytime operationswould be negligible. The additional operations during the daytime wouldconsist of computing the statistics from the initial state of requests setand providing them to the prediction model. The run of the second al-gorithm and the training of the prediction model might be done duringthe nighttime.This section presents the statistics computed from the benchmarkproblems and authors’ approach to creating a linear model, based onthose statistics, predicting the relative MEMSO and 2MPSO algorithmsperformance.
In order to create a set of features, allowing to discriminate betweendifferent benchmark problems, the authors have computed a set of statis-tics for each of the benchmarks. Those statistics may be divided into 2groups:requests set spatial features,requests set volume features.Within those groups, the following statistics have been computed:
Table 1.
Values of the input features computed for the initially known sets of requestsfrom Kilby’s and Montemanni’s benchmarks
Name µ x sd x skew x µ y sd y skew y µ s sd s skew s nc c50 0.53 0.32 -0.18 0.51 0.33 0.08 0.09 0.04 0.58 2.00c75 0.51 0.30 -0.08 0.41 0.29 0.27 0.12 0.06 0.34 4.00c100 0.48 0.27 0.35 0.46 0.27 0.53 0.08 0.05 0.71 1.00c100b 0.68 0.17 0.48 0.56 0.26 0.31 0.09 0.05 1.39 0.00c120 0.38 0.32 0.84 0.67 0.24 0.20 0.06 0.03 1.46 0.67c150 0.49 0.27 0.06 0.45 0.24 0.28 0.08 0.04 0.78 6.00c199 0.52 0.26 -0.04 0.47 0.24 0.19 0.09 0.04 0.52 9.00f71 0.44 0.23 -0.02 0.72 0.18 -0.17 0.04 0.05 1.87 0.00f134 0.61 0.29 -0.47 0.42 0.22 0.70 0.07 0.11 2.35 1.50tai75a 0.27 0.18 1.39 0.47 0.20 1.22 0.11 0.16 2.27 3.00tai75b 0.69 0.25 -1.14 0.42 0.22 -0.69 0.11 0.16 1.24 1.00tai75c 0.43 0.16 -0.93 0.49 0.25 -0.02 0.11 0.15 2.30 3.00tai75d 0.43 0.30 0.55 0.49 0.29 -0.23 0.09 0.13 1.55 0.25tai100a 0.39 0.30 0.61 0.57 0.25 -0.27 0.11 0.18 2.04 1.00tai100b 0.46 0.29 0.46 0.46 0.22 0.56 0.12 0.16 1.64 1.00tai100c 0.63 0.27 -0.91 0.50 0.20 -1.10 0.10 0.14 1.34 5.00tai100d 0.44 0.23 -0.39 0.46 0.21 0.06 0.10 0.17 2.37 2.00tai150a 0.47 0.29 0.17 0.55 0.31 -0.36 0.09 0.14 2.32 1.67tai150b 0.61 0.23 0.13 0.52 0.24 0.81 0.09 0.14 2.30 6.00tai150c 0.47 0.23 0.46 0.61 0.14 -0.06 0.09 0.13 1.80 6.00tai150d 0.47 0.34 -0.05 0.67 0.23 -0.88 0.12 0.16 1.41 0.14 Spatial features: – µ x , µ y - mean locations along coordinate system axes, – sd x , sd y - standard deviation of locations along coordinatesystem axes, – skew x , skew y - standardized skewness of locations along co-ordinate system axes, – k gap - gap statistic (estimated optimal number of clusters) forrequests locations,Volume features: – µ s - mean volume, – sd s - standard deviation of volume, – skew s - skewness of volume, – m v - minimum number of vehicles necessary to load all therequests. Particle Swarm Optimization hyper-heuristic for the Dynamic Vehicle Routing Problem In order to make the features comparable between different bench-marks, the spatial locations have been mapped to the [0 , × [0 ,
1] plain,requests volume has been divided by vehicles capacity and k gap has beencombined with the m v in the following way: nc = | − m v k gap | The larger values of nc suggest that the requests might not be easilydivided among the vehicles.The results of computing those features on the Kilby’s and Monte-manni’s benchmark set for the a priori available requests are presentedin Table 1. Listing 1.1.
Model trained on the full data set, with the variable selection usingAkaike Information Criterion, presented as an R output.
C o e f f i c i e n t s :E s t i m a t e Std . E r r o r t v a l u e Pr ( > | t | )( I n t e r c e p t ) 1 . 4 2 8 9 9 0 . 1 5 6 7 6 9 . 1 1 6 9 . 6 4 e − ∗∗∗ s d x 1 . 2 3 1 1 9 0 . 2 6 5 5 1 4 . 6 3 7 0 . 0 0 0 5 7 3 ∗∗∗ mu y − − − ∗∗∗ s d y − − ∗∗∗ skew y − − − − ∗∗∗ s k e w s 0 . 2 1 5 0 7 0 . 0 4 5 9 2 4 . 6 8 3 0 . 0 0 0 5 2 9 ∗∗∗−−− S i g n i f . c o d e s : 0 ∗∗∗ ∗∗ ∗ − s q u a r e d : 0 . 8 6 , A d j u s t e d R − s q u a r e d : 0 . 7 7 6 4F − s t a t i s t i c : 1 0 . 4 3 on 7 and 12 DF, p − v a l u e : 0 . 0 0 0 2 8 6In order to chose a proper optimization algorithm, the authors proposea simple linear regression model.To train such a model, for predicting the proper optimization algo-rithm for a given benchmark, the statistics presented in Tab. 1 has beenselected as an input variables and a ratio between the average MEMSOand the average 2MPSO result has been chosen as an output variable. Such approach is justified by the fact that it is more important to choosea proper algorithm when the difference between different algorithms per-formance is significant.Fitted linear model is used in a following way. The subset of thestatistics to be computed, described in Section 4.1, is identified withthe usage of a variable selection method (cf. Listing 1.1). When a newDVRP instance is considered, that subset of statistics, forming the inputvariables for the model, is computed over the requests of that instance.If the model predicts the relative result of MEMSO to 2MPSO to beless than 1, MEMSO is chosen as the optimization algorithm, otherwise2MPSO is chosen (see Fig. 2).
5. Results
Table 2.
Results of predicting the optimization algorithm from a leave-one-out cross-validation experiment. The table presents the minimum and the average valuesachieved by MEMSO and 2MPSO algorithms, marking the significantly better av-erage results with a gray background. The Hyper-heuristic column presents whichalgorithm has been chosen by a linear model, whether the choice has been appropri-ate and how much has been gained (or lost) with that choice of algorithm. For thebenchmark instances with significant difference between average results of MEMSOand 2MPSO the gain has been marked with a gray background and a loss with alight-gray background.
MEMSO T 4.15%c75
MEMSO T 1.80%tai75c
MEMSO T 1.81%tai75d 1445.58 1481.25
MEMSO T 2.23%tai100a 2211.30 2327.20
MEMSO T 4.87%tai100b
Particle Swarm Optimization hyper-heuristic for the Dynamic Vehicle Routing Problem To test the proposed approach a leave-one-out cross-validation exper-iment with a linear model predicting the relative performance of theMEMSO to 2MPSO has been performed. In each fold a linear modelhas been built on information from initial sets of requests from 20 out of21 benchmark problems. The relative algorithm performance has beencomputed as a ratio of the average results obtained by the MEMSO and2MPSO algorithms. Such approach simulated a performance of hyper-heuristic optimizing and training the model on some number of sub-sequently computed DVRP tasks. The average results were computedover 30 runs of MEMSO per benchmark and 20 runs of 2MPSO perbenchmark. The order of the total fitness function evaluations budgetfor each of the algorithms run has been equal to 10 . The details aboutparameter setting for MEMSO and 2MPSO experiments can be foundin [14] and [18], respectively.Additionally, in the model training phase some of the input variableshave been removed from the model in a step-wise mode, with the usage ofAkaike Information Criterion, in order to create a more general predictor.As an example, results of applying such procedure, to the model trainedon all 21 of the benchmarks, are presented as a listing of R output inListing 1.1. It can be observed (from the corresponding p -values) thatboth the spatial and the volume related features have been selected asinformative ones.The results from the cross-validation experiment are given in Table 2.The correct algorithm (the one with a better average result) has beenchosen in 14 out of 21 cases (all of them with statistically significant dif-ference between average algorithms results, verified by a t -test). Wrongpredictions, that have been made for the other 7 cases, have resulted ina significant loss of results quality only for 3 of them ( f71 , tai150a , and tai150b ). Therefore, the linear model achieved 82% accuracy for thesubset of 17 benchmarks with significant differences in the algorithmsaverage results (with 67% accuracy over the whole set of benchmarks).
6. Conclusions
Choosing an optimization algorithm for the DVRP on the basis ofthe characteristics of initial requests sets leads to the improvement ofthe results. Choosing between MEMSO and 2MPSO algorithms resultedin a 0.6% improvement in comparison with the 2MPSO average perfor-mance (with a maximum of 1.5% improvement if all the predictions wereaccurate) and 2.8% in comparison with the MEMSO performance. Thebest performance has been improved by 0.2% on average in comparisonwith 2MPSO and by 2.8% in comparison with MEMSO. The obtained results suggest that the proper benchmark features hasbeen selected and choosing an optimization algorithm for a given bench-mark problem is possible and may lead to results improvement. Both thespatial and the volume related features have been marked as significantin predicting relative performance.The future work should include an algorithm choice possibility duringthe optimization and not only at the beginning of the working day. Itmight be also beneficial to search for another set of characteristics al-lowing for a proper choice of the algorithm, in order to try eliminatingthe DVRP tasks for which the significantly worse algorithm has beenselected.
Acknowledgments
The research was financed by the National Science Centre in Polandgrant number DEC-2012/07/B/ST6/01527Project website: . eferences [1] E. Burke, G. Kendall, J. Newall, E. Hart, P. Ross, and S. Schulenburg. Hyper-heuristics: An emerging direction in modern search technology. Internationalseries in operations research and management science , pages 457–474, 2003.[2] N. Christofides and J. E. Beasley. The period routing problem.
Networks ,14(2):237–256, 1984.[3] C. W. Cleghorn and A. P. Engelbrecht. Particle swarm variants: standardizedconvergence analysis.
Swarm Intelligence , 9(2-3):177–203, 2015.[4] I. Cristian and Trelea. The particle swarm optimization algorithm: convergenceanalysis and parameter selection.
Information Processing Letters , 85(6):317–325,2003.[5] G. Croes. A method for solving traveling salesman problems.
Operations Res.6 , pages 791–812, 1958.[6] G. B. Dantzig and R. Ramser. The Truck Dispatching Problem.
ManagementScience , 6:80–91, 1959.[7] M. Elhassania, B. Jaouad, and E. A. Ahmed. A new hybrid algorithm to solvethe vehicle routing problem in the dynamic environment.
International Journalof Soft Computing , 8(5):327–334, 2013.[8] M. Elhassania, B. Jaouad, and E. A. Ahmed. Solving the dynamic vehicle rout-ing problem using genetic algorithms. In
Logistics and Operations Management(GOL), 2014 International Conference on , pages 62–69. IEEE, 2014.[9] M. L. Fisher and R. Jaikumar. A generalized assignment heuristic for vehiclerouting.
Networks , 11(2):109–124, 1981.[10] F. T. Hanshar and B. M. Ombuki-Berman. Dynamic vehicle routing usinggenetic algorithms.
Applied Intelligence , 27(1):89–99, Aug. 2007.[11] J. Kennedy and R. Eberhart. Particle Swarm Optimization.
Proceedings of IEEEInternational Conference on Neural Networks. IV , pages 1942–1948, 1995.[12] M. R. Khouadjia, E. Alba, L. Jourdan, and E.-G. Talbi. Multi-Swarm Op-timization for Dynamic Combinatorial Problems: A Case Study on DynamicVehicle Routing Problem. In
Swarm Intelligence , volume 6234 of
Lecture Notesin Computer Science , pages 227–238. Springer, Berlin / Heidelberg, 2010. [13] M. R. Khouadjia, B. Sarasola, E. Alba, L. Jourdan, and E.-G. Talbi. A com-parative study between dynamic adapted PSO and VNS for the vehicle routingproblem with dynamic requests.
Applied Soft Computing , 12(4):1426–1439, 2012.[14] M. R. Khouadjia, E.-G. Talbi, L. Jourdan, B. Sarasola, and E. Alba. Multi-environmental cooperative parallel metaheuristics for solving dynamic optimiza-tion problems.
Journal of Supercomputing , 63(3):836–853, 2013.[15] P. Kilby, P. Prosser, and P. Shaw. Dynamic VRPs: A Study of Scenarios, 1998.[16] R. Montemanni, L. Gambardella, A. Rizzoli, and A. Donati. A new algorithmfor a dynamic vehicle routing problem based on ant colony system.
Journal ofCombinatorial Optimization , 10:327–343, 2005.[17] M. Okulewicz and J. Ma´ndziuk. Application of Particle Swarm OptimizationAlgorithm to Dynamic Vehicle Routing Problem. In L. Rutkowski, M. Kory-tkowski, R. Scherer, R. Tadeusiewicz, L. Zadeh, and J. Zurada, editors,
ArtificialIntelligence and Soft Computing , volume 7895 of
Lecture Notes in Computer Sci-ence , pages 547–558. Springer Berlin Heidelberg, 2013.[18] M. Okulewicz and J. Ma´ndziuk. Two-Phase Multi-Swarm PSO and the Dy-namic Vehicle Routing Problem. In , pages 86–93, Orlando, Fl, USA, 2014.IEEE Press.[19] M. Okulewicz and J. Ma´ndziuk. Dynamic Vehicle Routing Problem: AMonte Carlo approach. In
Information Technologies: Research and TheirInterdisciplinary Applications 2015 , pages 119–138, Pu ltusk, Poland, 2015. http://phd.ipipan.waw.pl/pliki/mat_konferencyjne/12_ITRIA_2015_01.pdf .[20] V. Pillac, M. Gendreau, C. Gu´eret, and A. L. Medaglia. A review of dynamicvehicle routing problems.
European Journal of Operational Research , 225(1):1–11, 2013.[21] H. N. Psaraftis. Dynamic programming solution to the single vehicle many-to-many immediate request dial-a-ride problem.
Transportation Science , 14(2):130–154, 1980. cited By 198.[22] H. N. Psaraftis. Dynamic vehicle routing: Status and prospects.
Annals ofOperations Research , 61(1):143–164, 1995.[23] Y. Shi and R. Eberhart. A modified particle swarm optimizer.
Proceedingsof IEEE International Conference on Evolutionary Computation , pages 69–73,1998.[24] ´E. D. Taillard. Parallel iterative search methods for vehicle routing problems.