Mathematical Problems in Engineering

Mathematical Problems in Engineering / 2020 / Article

Research Article | Open Access

Volume 2020 |Article ID 7258780 | https://doi.org/10.1155/2020/7258780

H. Xue, "Adaptive Cultural Algorithm-Based Cuckoo Search for Time-Dependent Vehicle Routing Problem with Stochastic Customers Using Adaptive Fractional Kalman Speed Prediction", Mathematical Problems in Engineering, vol. 2020, Article ID 7258780, 18 pages, 2020. https://doi.org/10.1155/2020/7258780

Adaptive Cultural Algorithm-Based Cuckoo Search for Time-Dependent Vehicle Routing Problem with Stochastic Customers Using Adaptive Fractional Kalman Speed Prediction

Academic Editor: A. M. Bastos Pereira
Received24 Mar 2020
Revised27 Jun 2020
Accepted30 Jun 2020
Published24 Jul 2020

Abstract

For the Time-Dependent Vehicle Routing Problem with Stochastic Customers (TDVRPSC), an adaptive Cultural Algorithm-Based Cuckoo Search (CACS) has been proposed in this paper. The convergence of the new algorithm is proved. An adaptive fractional Kalman filter (AFKF) for traffic speed prediction is proposed. An adaptive mechanism for choosing the covariance of state noise is designed. Its mathematical process is proved. Several benchmark instances with different scales are tested, and new solutions are discovered, which are better than the published solutions. The effects of the parameters on the convergence and the results are studied. According to cargo weight of customers to be delivered, the customers can be divided into large, small, and retail customers. The algorithm is tested with fixed demand probability and also different customer types with stochastic demand. The traffic speeds in different business districts in Xiamen at different times are predicted by AFKF. The results show that AFKF has smaller prediction error and better prediction accuracy than fractional Kalman filter and Kalman filter. The effect of different fractional orders on prediction error is compared. The performance of the new algorithm is compared with that of the cultural algorithm and the Cuckoo Search. The result shows that the new algorithm can efficiently and effectively solve DTVRPSC and improve the accuracy of vehicle routing planning of time-varying actual urban traffic road.

1. Introduction

With the rapid development of the logistics industry, Vehicle Routing Problem (VRP) gets more and more attention. The traditional VRP uses static road network model and cannot accurately estimate the travelling time based on the actual changing traffic conditions. In the actual traffic applications, there are many uncertainties, such as random customer demand, random number of customers, random traffic speed, random conditions of road traffic, and rush hours. Therefore, Time-Dependent Vehicle Routing Problem with Stochastic Customers (TDVRPSC) is a stochastic optimization problem with extensive practical significance. It is important to study new intelligent algorithms to realize logistics automation under dynamic traffic information. The results of this paper can be applied to vehicle routing, artificial intelligence, and data filtering.

Guo studied dynamic multiobjective Vehicle Routing Problem with a corresponding carbon emission model and it was set as an optimization objective [1]. After finding optimal robust virtual routes for all customers by adopting multiobjective particle swarm optimization in the first phase, static vehicle routes for static customers are formed by removing all dynamic customers from robust virtual routes in next phase. The dynamically appearing customers append to be served according to their service time and the vehicles’ statues. Global vehicle routing optimization is triggered only when no suitable locations can be found for dynamic customers. A metric measuring the algorithms robustness is given. Guo proposed an adaptive immune clonal selection cultural algorithm [2]. Dual structure of cultural algorithm was adopted, and a hybrid selection strategy integrating selection and clonal selection was put forward. The proportion of population influenced by each selection method was adaptively adjusted according to implicit knowledge extracted from the evolution process. Guo proposed a novel multipopulation cultural algorithm adopting knowledge migration derived from human cultural interaction [3]. Knowledge extracted from the evolution process contained more effective information. By migrating knowledge among subpopulations at constant intervals, the algorithm realized more effective interaction with less communication cost. Geng presented a modified centralized algorithm based on particle swarm optimization (MCPSO) [4]. Three distributed algorithms were compared against three centralized algorithms. Results showed that MCPSO could be used as a benchmark for distributed algorithms for determining the performance gap.

Random VRP was first proposed by Stewart and Golden and has gotten more and more attention [5]. Bertsimas considered a natural probabilistic variation of the classical vehicle routing problem and used stochastic demands [6]. Bianchi et al. analyzed the performance of metaheuristics on the Vehicle Routing Problem with stochastic demands [7]. Bozorgi-Amiri et al. developed a multiple objective robust stochastic programming approach for disaster relief logistics under uncertainty [8]. Sungur et al. introduced a robust optimization approach to solve the VRP with demand uncertainty [9]. Tan et al. considered VRP with stochastic demand, where actual demand was revealed only when the vehicle arrived at the customer [10]. Ibarra-Rojas et al. studied a VRP to optimize accessibility based on six indicators [11]. Wollmer et al. used a cutting plane algorithm for solving the stochastic programming and proved the algorithm to be finite [12]. Kenyon and Morton used Monte Carlo sampling based solution procedure for stochastic VRP with large sample spaces [13]. The solution space of the classical NP-hard problem is of large scale and complex. Large-scale networks may be found in many fields and a crucial step is to estimate them from the data. Sparse reciprocal graphical models have been considered for learning large-scale networks using efficient algorithms. Alpago et al. proposed an identification procedure of a sparse graphical model associated with a Gaussian stationary stochastic process [14].

At present, perfect solutions can only be reached as close as possible within a reasonable run time. Therefore, new intelligent methods are required. In 2009, Yang and Deb proposed Cuckoo Search (CS) to solve optimization problems [15]. Mareli and Twala developed three Cuckoo Search algorithms based on dynamically changing switch parameters [16]. Ishak Boushaki et al. proposed a quantum chaotic Cuckoo Search algorithm for data clustering [17]. Ayoubi et al. proposed a nonhomogeneous Cuckoo Search algorithm for allocation of passive filters [18].

To accelerate the convergence speed, this paper adopts the double mechanism in cultural algorithm (CA) [19] and proposed a Cultural Algorithm Based Cuckoo Search (CACS). CACS can overcome the defects of slow convergence rate of CS, as well as the blind guidance of CA group space, which makes the whole population algorithm not easy to converge. Its adaptive mechanism can overcome the dependence of parameters on algorithm performance. In CACS, the individuals with high fitness form the corresponding culture of the group in the Belief space by accepting the rules. The representative individuals in the population, utilizing the outstanding individuals after adjustment, will further affect the evolution of next-generation population.

Kalman filter was proposed by Kalman in 1960 [20]. In recent years, it has been widely applied to adaptive control and system recognition [21]. With the presence of uncertainties, the robust Kalman filter based on the tau-divergence which accounts for uncertainties could be used. San Gultekin et al. considered the nonlinear Kalman filtering problem using α-divergence measures as optimization criteria [22]. In order to extend the application of Kalman filter from integer-order to fractional-order system, Solís-Pérez et al. presented a fractional-order Kalman filter according to the Riemann–Liouville definition [23]. Sun et al. used the extended Kalman filtering for fractional-order nonlinear discrete system with Lévy noises [24]. Liu et al. proposed a generalized fractional central difference Kalman filter for nonlinear discrete fractional dynamic systems [25].

In this paper, an adaptive fractional Kalman filter for traffic speed prediction is proposed to update and estimate the time-varying vehicle speed. An adaptive mechanism is adopted in computing the covariance of prediction error. Its mathematical process is proved. A solution of classic VRP cases that is better than those optimal solutions published before is found. Using AFKF and real-time urban traffic data, vehicle speed is predicted to adjust vehicle routing in time to avoid congested roads and improve the accuracy of VRP. The contributions of the paper are summarized as follows:(1)An adaptive Cultural Algorithm Based Cuckoo Search is proposed(2)An adaptive fractional Kalman filter for traffic speed prediction is proposed(3)The convergence of CACS and the mathematical process of AFKF are proved(4)New solutions of VRP are discovered, which are better than the published solutions(5)With real-time urban traffic data, vehicle speed is dynamically predicted to adjust the solution scheme of TDVRPSC in time

2. Mathematical Model of TDVRPSC

2.1. Mathematical Model

Stochastic time-dependent road network can be represented by a directed network . represents a set of customer nodes, and represents a set of road segments. . represents the travel speed from customer node i to j during time period .

The goal is to obtain the shortest expected travel time and service time. denotes the customer demand probability. The total number of vehicles is K. denotes the travelling time from customer i to j. denotes whether vehicle l is responsible for the route from customer j to i. denotes whether customer i is serviced by vehicle l. denotes the service time of customer i by vehicle l. is the demand of the customer i. is the total load of the vehicle l. The model can be described as follows:

Hypothesis is stated as follows:(1)Each customer is serviced exactly once. The loading and unloading services of each customer are completed only once. Equations (2) and (3) ensure that each customer is served exactly once.(2)Each customer is loaded immediately after being visited by a vehicle. The total time includes travelling time and service time with no other times wasted.(3)The total weight of goods demanded by customers cannot exceed the vehicle load servicing them. Equation (6) ensures the total weight of goods demanded by customers cannot exceed the vehicle load servicing them.(4)Each customer is serviced by a vehicle only if he is visited on the route of the vehicle. Equation (5) ensures a customer is serviced by a vehicle only if he is visited on the route of the vehicle.(5)Each customer can only be serviced by exactly one vehicle. Equation (4) ensures each customer is serviced exactly one vehicle.(6)The positions of all customers are known.(7)The position of the depot is known. The source of vehicles is determinant.(8)The distances from the depot to each customer are known.(9)The customer demand probability is a stochastic real number between zero and one.(10)The travelling time between two customers is stochastic time-dependent due to the time-varying speed and real-time traffic conditions on roads.

The existence of at least one solution of the above problem is proved in literature [12, 13]. A cutting plane algorithm for solving the above stochastic programming was proved to be finite [8]. The above stochastic objective function (1) was proved to be a convex function on the convex hull of the set of decision vectors satisfying constraints (2)–(8) [9]. Since the local minimum of convex optimization problem is the global minimum, the local minimum solution of the above problem will also be the global minimum.

2.2. Adaptive Fractional-Order Kalman Filter for Speed Prediction

The distance between customer and is represented by . denotes the sampling time interval. The calculation procedure for travelling time is shown in Figure 1. is updated and estimated by Kalman filter algorithm.

The fractional-order difference is formulated aswhere is the fractional order. is the sampling interval. is the number of samples. The fractional-order model of speed prediction is given bywhere is the state variable of the traffic speed. is the measurement output. is the system state noise. is the measurement noise. is calculated aswhere is the measurement sequence containing measurement output . denotes the covariance of prediction error at time instant , which is calculated as

denotes the covariance of system noise at the time instant , which is calculated as

denotes the covariance of measurement noise at the time instant , which is calculated as

denotes the state estimation at time instant , which is calculated as

is the covariance of estimation error at the time instant , which is calculated as

The covariance of state noise is designed as follows:where is a gain coefficient.

When the covariance of state noise is small, the convergence speed is slow and the convergence accuracy is high. If the state changes too fast in a short time, the slow convergence rate will easily lead to tracking failure. When the covariance of state noise is large, the convergence speed is fast and the accuracy is low. Before reaching a stable state, choosing a larger covariance of state noise can accelerate the convergence speed. After reaching a stable state, the convergence accuracy can be improved by choosing smaller covariance of state noise.

In standard Kalman filter, is a constant. However, it is difficult to determine an optimal gain in the whole process. When the covariance of state noise is small, the number of estimations will update slowly if is too small. When the error changes sharply, the estimation will oscillate and even the convergence process will be unstable if is too large. To accelerate the convergence process, an adaptive gain is used to make the estimation process adjusted reasonably and adapt to different data characteristics. increases when the covariance of state noise is small and decreases when the covariance of state noise is large. According to the noise covariance information, is calculated dynamically. Therefore, an adaptive mechanism for choosing is designed as follows:

The Kalman gain is calculated as

The fractional Kalman filter is formulated as

Assumption 1. . is not related to .

Assumption 2. . is not related to and .

Theorem 1. For the fractional-order system with system state noise and measurement noise described by formulas (11)–(13), AFKF is given by formulas (23)–(27).

Proof. Substituting (11) and (12) into (14), one can obtainFrom Assumption 1, one can obtainSubstituting (29) into (28), one can obtainSubstituting (18) into (29), one can obtainSubstituting (18) into (31), one can obtainFormula (23) is obtained.
Substituting (11) and (12) into (15), one can obtainSubstituting (32) into (33), one can obtainSubstituting (19) and (20) into (34), one can obtainFormula (27) is obtained.
The cost function is defined byDifferentiating (36) with respect to and making the result equal to zero, one can obtainFrom (37), one can obtainSubstituting (22) into (38) yieldsFormula (25) is obtained.
Substituting (39) into (19) yieldsSubstituting (13) into (40) yieldsSubstituting (15) and (17) into (41) yieldsFrom Assumption 2, one can obtainSubstituting (43) into (42) yieldsFormula (26) is obtained.
Substituting (13) into (17) yieldsSubstituting (15) into (45) yieldsFrom Assumption 2, one can obtainSubstituting (15) into (45) yieldsFormula (24) is obtained.

3.1. Encoding Method

Construct the mapping relationship from the individual to the solution decision variable . is composed of randomly generated (N + K − 1)-dimensional real numbers. Take Table 1 for example. Suppose N is 4 and K is 3. N + K − 1 = 6.


j

1−0.992
21.805
33.014
4−0.721
5−1.206
62.163

is the largest dimension, so , is the second largest dimension, so , and so on. Finally, the whole mapped individual is .

Then, are regarded as the boundary lines of different vehicles. Customers within the same boundary lines are assigned to the same vehicle in turn. So and are regarded as boundary lines. So the mapped solution decision variable can be written as . This means that the second customer is assigned to the first vehicle, the fourth and the first customers are assigned to the second vehicle, and the third customer is assigned to the third vehicle.

3.2. Cuckoo Search

In CS, the position value of each cuckoo nest corresponds to each solution to VRP. Its goal is to search for the best cuckoo nest.

Lévy flight is a random walk strategy, whose step length can be determined by the Lévy distribution. When the new is discovered to be solution-suitable, a Lévy flight can be used as follows:where i represents the ith cuckoo. represents the step size related to the scales of the problem. is the index coefficient. An adaptive method can also be used:where is a constant. The product means entry-wise multiplications. The step length of a random walk by the Lévy flight is computed from a Lévy distribution:

The steps essentially form a random walk process with a power-law step-length distribution with a heavy tail. Some new solutions should be generated by Lévy walk around the best solution obtained so far to speed up the local search. The locations of a new solution may be far enough from the current best solution; this will make sure the system will not be trapped in a local optimum. When the length of random step is determined by the maximum step length in the Lévy distribution, the Lévy flight can effectively provide a random walk model.

To avoid parameters dependence, an adaptive mechanism is introduced. The step size of Lévy flight and discovery probability are computed as follows:where denotes the maximum step size of Lévy flight and denotes the minimum step size of Lévy flight. denotes the maximum fitness of all nests, and denotes the minimum fitness of all nests. denotes the maximum discovery probability, and denotes the minimum discovery probability. denotes the current iteration, and denotes the maximum iterations number. As the fitness of a nest decreases, its step size of Lévy flight will increase. As the time of iteration increases, the discovery probability of nests will decrease.

3.3. Cultural Algorithm

CA is a dual evolution mechanism proposed by Reynolds in 1994. For VRP, situation knowledge represents the optimization goal, and normative knowledge represents constraint conditions. The individuals gained in the evolution process are passed to Belief space by Accept function. Then the Belief space is updated through Update function. The Influence function uses the useful information of Belief space to guide the evolution process of population space. The main structure of the cultural algorithm is shown in Figure 2.

In the Belief space, denotes the upper bound of the individual, and denotes the lower bound of the individual. is lower limit of the fitness value of , and is upper limit of the fitness value of . The updating rules of situation knowledge are as follows:where denotes the optimal individual in generation and denotes the fitness function. The updating rules of standard knowledge are

3.4. CACS Process

To accelerate the convergence speed of CS, this paper adopts the double mechanism of the CA. The individuals with high fitness form the culture in the Belief space through acceptance rules. The representative individuals further influence the evolution of next population through the outstanding individuals adjusted. The process of CACS with vehicle speed estimation is shown in Figure 3.

In this paper, the adaptive factor strategy is introduced. Through the dynamic adjustment of the factor to the optimization step size, the algorithm can search with adaptive step size in the early stage of evolution and improve the convergence speed of the algorithm. In the late stage of evolution, it can quickly tend to the optimal solution and further improve the convergence speed of the algorithm. This reduces the computational cost of the algorithm.

3.5. Convergence of CACS

The solution generated in the kth iteration is denoted as . The optimal group state set is denoted as H. The optimal state set is denoted as R. All the states of all the cuckoos form a state space for the group denoted as Q.

Lemma 1. (see [26]). Assume that the objective function f is a measurable function. The feasible solution space A is a measurable subset on . is a set of global optimal solutions. denotes the probability that the solution belongs to the optimal solutions set . For any Borei subset B, its Lebesgue measure . is the probability of obtaining B by measuring .. Thus,

Lemma 2. The finite homogeneous Markov chain starting from any nonrecurrent states will reach its recurrent state with probability of 1.

Lemma 3. Suppose Markov chain has a nonempty set C and there is not any nonempty closed set D satisfying . Then when , .

Theorem 2. The Markov chain represented by CS is finite and homogeneous.

Proof. The search space of the optimization algorithm is limited, and the cuckoo nest state is limited. So the state space of nest position is limited. Since the population state of nest positions is composed of m nest positions and m is a finite positive integer, the population state of nest position is also limited.
From the process of CS, the transfer probability of any nest position is only related to the state of the previous state . denotes the transfer operation. Therefore, the population states of nest position have Markov property.
The transfer probability of any nest position is not related to the time instant k. That is to say, the population states sequence of nest location is homogeneous. Therefore, the population states sequence of cuckoo nest location is a finite homogeneous Markov chain.

Theorem 3. When the number of iterations approaches becoming sufficiently large, the group state sequence will converge to the optimal state solution set H.

Proof. , . From the process of CS, it is obtained that . For the state sequence {y (k); k > 0}, the state set R of optimal nest position is a closed set on the state space Y.
Since R is a closed set on Y, according to the definition of closed set, the state set H of the optimal nest position group is a closed set on the group state space Q.
Assume that there exists a closed set B in the group space Q, such that . , . Then, , . Thus, it is obtained that . So B is not a closed set, which contradicts the hypothesis. Therefore, there is no nonempty closed set outside H in the nest position state space Q such that .
From Lemma 3, when the number of iterations approaches becoming sufficiently large, the group state sequence will converge to the optimal state solution set H.

Each iteration of CS algorithm preserves the optimal position of the population and guarantees the nonincremental fitness. According to Theorem 3, the population sequence of cuckoo nest position will enter the optimal state after continuous and infinite iterations. Therefore, the possibility that the global optimal solution cannot be searched continuously and infinitely is 0. So the CS algorithm converges to the global optimum.

Theorem 4. The random variation of the state in Belief space of CACS belongs to finite homogeneous Markov chain.

Proof. In CACS, the Belief space is used to record the evolutionary knowledge in the search process, and CACS searches for the optimal value in discrete and limited space. So the Belief space used for recording the evolutionary knowledge is also limited. Meanwhile, the state transition probability in the Belief space is only related to its current state and has nothing to do with time instant k. Therefore, the random variation of the state in the Belief space of CACS belongs to the finite homogeneous Markov chain.

Theorem 5. CACS has global convergence with probability of 1.

Proof. Denote the global knowledge state space in the Belief space as . Denote the local knowledge state space as . The whole state space of Belief space is , in which denotes the global knowledge set and denotes the local set of population space.
is divided into two subspaces, namely, recurrent state space and nonrecurrent state space . The recurrent state space can be expressed as follows:and thus and .
If there is a knowledge state , according to the definition of global knowledge, it will not transfer to . So it is a closed set. The state in is interlinked. So the state in is a recurrent state. If , each individual will move throughout the whole solution space due to evolution. So the convergence probability of the evolutionary population to the global optimal solution is greater than 0. That means the transition probability of from to is greater than 0. Therefore, the state in is a nonrecurrent state. According to Lemma 2, the evolution of knowledge in Belief space from any nonrecurrent state will always transfer to the recurrent state with probability of 1, that is, converge to the global optimal solution according to probability of 1.

4. Example and Analysis

4.1. Experimental Results with Deterministic Demand Probability

The CACS algorithm is tested on Intel® Core™ i3-4150T CPU @3.00 GHz, 64-bit operating system with memory of 4.00 GB, and x64-based processor. The VRP library used is published on the Networking and Emerging Optimization homepage [27]. The capacitated VRP instance library was selected. The first layer index is the first letter of word in the instance library. The second layer index is the number of customers, and the third layer index is the number of vehicles. Take A, B, P, and E sets as examples. Set A has 27 VRP instances from A-n32-k5 to A-n80-k10 proposed by Augerat et al. Set B has 23 VRP instances from B-n31-k5 to B-n78-k10 proposed by Augerat et al. Set P has 24 VRP instances from P-n16-k8 to P-n101-k4 proposed by Augerat et al. Set E has 15 VRP instances from E-n13-k4 to E-n101-k14 proposed by Augerat et al.

Consider the case when the arrival probability is 1. Table 2 gives the shortest path length, average iteration times, and average CPU computing time of each instance.


InstanceA-n32-k5A-n39-k6B-n31-k5P-n16-k8P-n19-k2

N3239311619
K56582
Optimal solution (m)784.612831.2696662.3275451.9471212.6569
Average iteration number427648329193308
Average CPU computing time (s)18.968827.15636.81253.29696.7188
Published optimal solution (m)784831672450212
Improvement ratio (%)001.4393600

For instance, B-n31-k5 algorithm has found a better solution than the published optimal solution, whose scheme is shown in Table 3.


RoutePublished optimal solutionNew solution
Travelling sequenceDistanceCapacityTravelling sequenceDistanceCapacity

130 ⟶ 23 ⟶ 8 ⟶ 12 ⟶ 28 ⟶ 26109.6811977 ⟶ 30 ⟶ 5 ⟶ 25 ⟶ 18 ⟶ 16 ⟶ 21119.237289397
221 ⟶ 16 ⟶ 18 ⟶ 25 ⟶ 5 ⟶ 4 ⟶ 29112.99158612 ⟶ 24 ⟶ 19 ⟶ 1 ⟶ 3 ⟶ 6 ⟶ 9102.365980299
37 ⟶ 17 ⟶ 13 ⟶ 6 ⟶ 9 ⟶ 22110.82879620 ⟶ 27 ⟶ 10 ⟶ 2 ⟶ 14 ⟶ 15 ⟶ 11 ⟶ 28 ⟶ 8240.508851895
420 ⟶ 27 ⟶ 10 ⟶ 2183.5752384 ⟶ 29 ⟶ 22 ⟶ 13 ⟶ 17 ⟶ 23108.715642299
514 ⟶ 15 ⟶ 11 ⟶ 24 ⟶ 19 ⟶ 1 ⟶ 3155.4298952698.0612053822
Total672.506412668.889412

Figure 4 shows the shortest path of 70 nodes when demand probability is 1. The horizontal axis represents the horizontal axis of a city map with the unit of meter. The vertical axis represents the longitudinal axis of a city map with the unit of meter.

Figure 5 shows the objective function convergence curves of 70 nodes when demand probability is 1. The horizontal axis represents iteration times. The vertical axis represents the best individual in each iteration.

4.2. Parameter Influence Analysis

Figure 6 shows the average results of E-n22-k4, B-n31-k5, and B-n64-k9 after 20 times calculation under different demand probability with IT 100. The updating rate of situational knowledge is 0.3. The updating rate of standard knowledge is 0.5. The population size is 50. is 1.

Figure 7 shows average results of 20 times calculation under different with IT 100. is 0.45. is 0.5 and is 50.

Figure 8 shows average results of 20 times calculation under different with IT 100. is 0.45. is 0.45 and is 50.

Figure 9 shows average results of 20 times calculation under different with IT 100. is 0.45. is 0.3 and is 0.5.

From the above figures, when m is larger, the optimal solution is better. Increasing m is advantageous to searching for better solutions. However, when m reaches a certain value, the optimal solution is almost unchanged as m increases. In this instances, the optimal parameters are , , , and .

The sensitivity analysis of iterations times is shown in Table 4. In the case of E-n22-k4, when IT is too small, the calculation results still need to be improved and is large. When IT reaches 130, remains close to zero as IT increases. The results are almost stable.


Nmax100110120130140

Mean379.3766377.4962375.8675375.1023375.1023
Std.0.76250.21660.10260.01250.0063
ΔNmax1010101010
ΔMean1.88041.62870.765200
ΔMean/ΔNmax0.188040.112870.0765200

We have also tried to vary the number of host nests and the discovery probability in Cuckoo Search algorithm. We have used n = 5, 10, 15, 20, 50, 100, 150, 250, and 500 and  = 0.01, 0.05, 0.1, 0.15, 0.2, 0.25, 0.4, and 0.5. We found that n = 15 and  = 0.25 are sufficient for most optimization problems. Results and analysis also imply that the convergence rate of Cuckoo Search algorithm is not sensitive to the parameters used.

Apart from the number of host nests, there is one parameter . This also means that we do not have to fine-tune these parameters for a specific problem. Subsequently, Cuckoo Search algorithm is more generic and robust for many optimization problems. There are fewer parameters to be fine-tuned in Cuckoo Search than in other metaheuristic algorithms.

4.3. VRPSC Experiments with Fixed Demand Probability

is set as different real numbers from 0.1 to 0.9. The parameters are set as follows: , , , and . The initial vehicle speeds are set as stochastic real numbers within the range of [1, 6]. Table 5 gives the average results of shortest expected total travelling time and servicing time for P-n16-k8 case with different IT and different under 20 times calculation.


100200300400500

0.147.01664747.01664747.01664747.01664747.016647
0.265.35669665.35668065.35668065.35668065.356681
0.374.81729074.81697874.81697574.81697774.816984
0.481.35662781.35402281.35401081.35401781.354064
0.586.76676686.75332886.75322086.75316486.753328
0.691.63513591.58497891.58353291.58332291.583318
0.796.18311196.04008496.02895496.02695096.026981
0.8100.469282100.163211100.110735100.106147100.104906
0.9104.306277103.887785103.744337103.739400103.737179

Table 5 shows that, with the increasing of , the objective function also increases. With the increasing of IT, the objective function reduces. However, when IT reaches 300, even if IT continues to increase, the shortest expected total travelling time and servicing time will not be reduced and will remain almost the same.

4.4. VRPSC Experiments with Distinguishing Customer Characteristics

Customers can be divided into large customers, small and medium customers, and retail customers according to the weight of the goods they order. Consider the influence of customer property in these three cases. In the case of B-n31-k5, the customers who order goods weighing in the range of [19, 25] belong to large customers. The customers who order goods weighing in the range of [9, 18] belong to small and medium customers. The customers who order goods weighing in the range of [2, 8] belong to small and retail customers. The parameters are set as follows: , , , , and . The initial vehicle speeds are set as stochastic real numbers within the range of [1, 6].

Table 6 shows the optimal solution, average solution, and standard deviation in three cases of different customer types after 20 times calculation, when the demand probability is 0.8 or 0.9. In Table 6, 1st represents the result for the first calculation, 2nd denotes the results for the second calculation, and so on.


Probability0.90.8
CustomerLargeSmallRetailLargeSmallRetail

1st596.155638.2462656.2919585.8801590.4435638.0412
2nd594.4813637.9749651.5232585.6425594.5347636.5034
3rd590.7583635.9657653.6638585.8801594.5347637.8991
4th596.155638.7321659.7856587.4265594.5347636.5034
5th590.7583639.7165653.4783582.6757590.4435634.8529
6th590.7583637.8215655.5595588.2996590.4435636.5034
7th590.7583637.9749657.6003582.6757590.4435634.8529
8th596.155637.8215656.7167588.2996594.5347636.8692
9th590.7583637.9749651.375582.6757590.4435636.5034
10th596.155638.6415654.8217588.2996594.5347636.5034
Ave.593.2893638.087655.0816585.7755592.4891636.5032
Std.2.7141530.9476772.6643242.3735792.1562521.04854

Table 6 shows when the large customers have least objective function. The small and medium customers have second shortest objective function. The retail customers have the longest objective function. This is because when the vehicle does not need to serve large customers, it can freely and flexibly arrange its routes. When vehicles need to serve large customers, with the limit of maximum carrying capacity, the choice of driving routes has more constraints. With the increase of , the objective function also increases. The objective function with of 0.9 is greater than that with of 0.8.

4.5. Comparative Analysis

Table 6 gives performance of CACS, CA, and CS with different . Since CACS uses CA and CS, it is more complex than these two algorithms, while its time cost and storage cost are slightly higher compared to the latter two. However, the rapid development of computer hardware technology makes this difference for solving problems of general scale not apparent. From the experimental results, CACS performs better than CA, CS, and ACO.

To avoid sensitivity to parameters selection, which is one of the main limitations in swarm intelligence algorithms, adaptive mechanism is adopted into CACS. ACO [28] used for comparison also proposed machine learning strategies to control the parameter adaptation. The objective function by different and different algorithms is listed in Table 7.


CACSCACSACO
ResultsRate (%)ResultsRate (%)ResultsRate (%)

0.147.01664747.1056540.1889547.1341360.2492747.0910610.15802
0.265.35668165.3745870.0273965.3943740.0576465.3647340.01232
0.374.81698474.8348670.0239074.9212100.1391174.9014230.11273
0.481.35406481.5048640.1850281.4314570.0950481.4017630.05860
0.586.75332886.8912110.1586886.7931340.0458686.8014750.05547
0.691.58331891.7163010.1449991.6137060.0331791.6546010.07777
0.796.02698196.1053460.0815496.0713410.0461796.0578200.03210
0.8100.104906100.1159630.01104100.1100360.00512100.1134670.00855
0.9103.737179103.9213470.17722103.8641360.12223103.9163160.12223

4.6. TDVRPSC Using AFKF Speed Prediction

Based on the technology of AFKF, the highway data analysis platform is established. After extracting the information and constructing appropriate traffic prediction model, the future vehicle speed in different business districts can be effectively estimated. When the speed of a road is detected to be below the threshold, it can be inferred that this road is congested due to accidents or other reasons. Rapid identification of congested roads is completed. This is also useful for emergencies. When traffic congestion is detected in a certain road, decision support and vehicle dispatching optimization can adjust VRP scheme in time and avoid congested roads, so as to reduce transportation time and cost. In the field of real-time traffic forecasting, the fast information processing and mining ability of large data have very high reliability for real-time VRP scheduling. The accuracy of VRP scheduling is improved for the time-varying actual urban traffic road.

Table 8 lists the traffic speed and congestion delay index in different business districts of Xiamen at 16:25 on February 17, 2019 [29]. Congestion delay index is used to evaluate the degree of urban congestion. It is the ratio of the actual travel time to the free travel time of a trip.


RankingBusiness districtCongestion delay indexSpeed (km/h)

1Siming1.5827.21
2Huli1.3336.46
3Haicang1.2544.82
4Jimei1.2441.80
5Tongan1.1347.74

Table 9 lists the traffic speed in different business districts of Xiamen from 16:30 to 20:00 on February 17, 2019. The unit of speed is km/h.


TimeSimingHuliHaicangJimeiTongan

16:3027.2536.8845.1342.0346.89
16:4027.3336.3842.0146.1946.78
16:5027.8536.2545.3942.3446.05
17:0028.0736.5345.7242.5245.55
17:1028.8336.6942.3845.9945.89
17:2028.4636.4641.7045.7844.90
17:3028.2736.7140.9342.9245.64
17:4027.5736.6340.9242.8244.83
17:5026.8841.8136.3341.1843.01
18:0026.4140.5636.0942.4641.58
18:1026.6840.1335.6842.1441.38
18:2026.8140.2535.4140.8444.56
18:3026.5935.3741.6140.7645.95
18:4027.0441.4336.1441.1346.18
18:5027.5436.5243.2941.6646.65
19:0028.4337.0645.6042.4747.07
19:1029.0737.0746.2343.4747.92
19:2029.9437.7846.9343.547.93
19:3030.2138.0646.6743.8348.18
19:4031.0438.2147.0144.2248.78
19:5031.2238.4647.0243.9548.50
20:0031.6938.9447.4044.0249.10

Figure 10 shows Xiamen city map at 16:25 on February 17, 2019. Different business districts are marked with their congestion delay indexes.

Several branches of Yuantong Express in Xiamen city distributed in different business districts are used. Table 9 lists their business district, latitude, longitude, and the Mercator coordination converted. The No. 0 branch is the head office. VRP task is to distribute goods to all other branches.

4.7. Accuracy of Prediction Model

Figure 11 shows the traffic speed prediction in Huli district using AFKF prediction from 16:30 to 20:00 on February 17, 2019. The horizontal axis represents time. The vertical axis represents the traffic speed in km/h. The star marking connected with solid wire represents the real speed. The joint spider connected with dotted lines represents the predicted speed.

Table 10 shows the prediction mean square error (MSE) and mean absolute percentage error (MAPE) by AFKF, fractional Kalman filter (FKF), and Kalman filter (KF). Table 11 shows AFKF has smaller prediction error and better prediction accuracy than FKF and KF.


Business districtAFKFFKFKF
MSEMAPEMSEMAPEMSEMAPE

Siming1.44433.87731.98884.59952.62755.2535
Huli2.76323.49292.78233.61893.10273.6632
Haicang2.29123.42452.33343.43222.56533.4551
Jimei1.75562.29081.76632.49312.66642.736
Tongan1.49512.07321.94252.21981.92842.4390


No.AreaLongitude (degree)Latitude (degree)Vertical coordinate (m)Horizontal ordinate (m)

0Siming118.07654824.449664006147762739976
1Siming118.08607424.483964006163052745670
2Siming118.13810824.477468006252982744930
3Siming118.1701624.482423006298182745249
4Huli118.14282224.514657006254372751505
5Huli118.17888724.504187006322212749519
6Huli118.14729724.488977006267372747251
7Huli118.14467624.521974006259592752476
8Jimei118.10840724.610044006201222768477
9Jimei118.1046424.608536006190482769188
10Tongan118.15615524.726039006277082790657
11Tongan118.14663524.697452006261892785614
12Tongan118.13816224.708646006249352787738
13Haicang118.01761924.531057006048642753869
14Haicang118.04275824.498528006086092748789
15Haicang118.01589924.530205006043892753588

Figure 12 shows that the prediction errors vary with different fractional orders of Kalman filter. The horizontal axis represents fractional order. The vertical axis represents the prediction error of traffic speed in km/h. The five-pointed star represents the mean relative error. The asterisk represents the mean square error. For certain business districts, there exists a most suitable fractional order for reducing its prediction error. Compared with integer-order filters, fractional-order Kalman filters have wider parameter selection range and higher prediction accuracy.

If the prediction model is not accurate, wrong route selection may be made. For example, the traffic flow is busy in Siming district and its vehicle speed is relatively slow. If the vehicle selects Siming district with no accurate speed prediction, it will take more travelling time, even if the distance is short.

4.8. Vehicle Routing Solution

With four vehicles, TDVRPSC schemes including the route number, customer visiting sequence, total route length, travelling starting time, travelling ending time, and total travelling time are shown in Table 12. The parameter settings are the same as those in Section 4.5.


RouteCustomer visiting sequenceLength (km)StartEndTime (min)

10 ⟶ 8 ⟶ 11 ⟶ 10 ⟶ 12 ⟶ 9 ⟶ 14 ⟶ 0109.603316:3018:40130
20 ⟶ 6 ⟶ 4 ⟶ 7 ⟶ 5 ⟶ 3 ⟶ 047.314616:3017:5787
30 ⟶ 2 ⟶ 1 ⟶ 026.549016:3017:2656
40 ⟶ 13 ⟶ 15 ⟶ 034.740716:3017:2757

Table 12 gives the minimal expected total travelling time with different , different departure time, and different algorithms. Its performance is compared with CA, CS, and ACO. The parameter settings of the algorithms are the same as those in Section 4.5.

Table 13 shows CACS is better than CA, CS, and ACO algorithms in solving TDVRPSC. With the increase of , the minimal total expected time also increases. The total time also depends on the departure time. When the departure time is not the rush hour with traffic congestion, the total time is lower. When the departure time is the rush hour with traffic congestion, the minimal total expected time will increase.


(%)Departure timeCACSCACSACO

2516:30223228225226
17:30231237234236

5016:30276282279280
17:30286281277279

7516:30312318314317
17:30324329326328

10016:30330337333336
17:30343349345347

Compared with other kinds of cultural algorithms, such as adaptive immune clonal selection cultural algorithm, the convergence rate of Cuckoo Search algorithm is not sensitive to the parameters used. Therefore, adaptive Cultural Algorithm Based Cuckoo Search has fewer parameters to be fine-tuned in Cuckoo Search and is more robust for optimization problems.

5. Conclusion

The main contributions of this paper are summarized as follows:(1)An adaptive CACS has been proposed and its convergence is proved.(2)AFKF for traffic speed prediction is proposed. An adaptive mechanism is adopted with the covariance of prediction error. Its mathematical process is proved.(3)New solutions of VRP are discovered, which are better than the published solutions.(4)CACS and AFKF are used for TDVRPSC with real-time urban vehicle data. Traffic speed is predicted and vehicle routing can be adjusted in time to avoid congested roads to reduce transportation time and cost.

The next step is to study better TDVRPSC models and design new methods to deal with uncertainties.

Data Availability

The data used to support the findings of this study are included within the article.

Conflicts of Interest

The author declares that there are no conflicts of interest.

Acknowledgments

This work was supported in part by the National Natural Science Foundation of China under Grant 51579114, the Natural Science Foundation of Fujian Province under Grant 2018J05085, and High Level Scientific Research Project of Transportation Engineering (202003).

References

  1. Y.-N. Guo, J. Cheng, S. Luo, D. Gong, and Y. Xue, “Robust dynamic multi-objective vehicle routing optimization method,” IEEE/ACM Transactions on Computational Biology and Bioinformatics, vol. 15, no. 6, pp. 1891–1903, 2018. View at: Publisher Site | Google Scholar
  2. Y.-n. Guo, H. Wang, and J. Cheng, “Adaptive immune clonal selection cultural algorithm,” Acta Electronica Sinica, vol. 38, no. 4, pp. 966–972, 2010. View at: Google Scholar
  3. Y.-n. Guo, J. Cheng, Y.-y. Cao, and Y. Lin, “A novel multi-population cultural algorithm adopting knowledge migration,” Soft Computing, vol. 15, no. 5, pp. 897–905, 2011. View at: Publisher Site | Google Scholar
  4. N. Geng, Q. Meng, D. Gong, P. W. H. Chung, and H. Chung, “How good are distributed allocation algorithms for solving urban search and rescue problems? a comparative study with centralized algorithms,” IEEE Transactions on Automation Science and Engineering, vol. 16, no. 1, pp. 478–485, 2019. View at: Publisher Site | Google Scholar
  5. W. R. Stewart and B. L. Golden, “Stochastic vehicle routing: a comprehensive approach,” European Journal of Operational Research, vol. 14, no. 4, pp. 371–385, 1983. View at: Publisher Site | Google Scholar
  6. D. J. Bertsimas, “A vehicle routing problem with stochastic demand,” Operations Research, vol. 40, no. 3, pp. 574–585, 1992. View at: Publisher Site | Google Scholar
  7. L. Bianchi, M. Birattari, M. Chiarandini et al., “Hybrid metaheuristics for the vehicle routing problem with stochastic demands,” Journal of Mathematical Modelling and Algorithms, vol. 5, no. 1, pp. 91–110, 2006. View at: Publisher Site | Google Scholar
  8. A. Bozorgi-Amiri, M. S. Jabalameli, and S. M. J. Mirzapour Al-e-Hashem, “A multi-objective robust stochastic programming model for disaster relief logistics under uncertainty,” OR Spectrum, vol. 35, no. 4, pp. 905–933, 2013. View at: Publisher Site | Google Scholar
  9. I. Sungur, F. Ordóñez, and M. Dessouky, “A robust optimization approach for the capacitated vehicle routing problem with demand uncertainty,” IIE Transactions, vol. 40, no. 5, pp. 509–523, 2008. View at: Publisher Site | Google Scholar
  10. K. C. Tan, C. Y. Cheong, and C. K. Goh, “Solving multiobjective vehicle routing problem with stochastic demand via evolutionary computation,” European Journal of Operational Research, vol. 177, no. 2, pp. 813–839, 2007. View at: Publisher Site | Google Scholar
  11. O. J. Ibarra-Rojas, L. Hernandez, and L. Ozuna, “The accessibility vehicle routing problem,” Journal of Cleaner Production, vol. 172, pp. 1514–1528, 2018. View at: Publisher Site | Google Scholar
  12. R. D. Wollmer, “Two stage linear programming under uncertainty with 0-1 integer first stage variables,” Mathematical Programming, vol. 19, no. 1, pp. 279–288, 1980. View at: Publisher Site | Google Scholar
  13. A. S. Kenyon and D. P. Morton, “Stochastic vehicle routing with random travel times,” Transportation Science, vol. 37, no. 1, pp. 69–82, 2003. View at: Publisher Site | Google Scholar
  14. D. Alpago, M. Zorzi, and A. Ferrante, “Identification of sparse reciprocal graphical models,” IEEE Control Systems Letters, vol. 2, no. 4, pp. 2475–1456, 2018. View at: Publisher Site | Google Scholar
  15. X. S. Yang and S. Deb, “Cuckoo search via Lévy flights,” World Congress on Nature and Biologically Inspired Computing (NaBIC), IEEE Publications, New York, NY, USA, 2009. View at: Google Scholar
  16. M. Mareli and B. Twala, “An adaptive Cuckoo search algorithm for optimisation,” Applied Computing and Informatics, vol. 14, no. 2, pp. 107–115, 2018. View at: Publisher Site | Google Scholar
  17. S. Ishak Boushaki, N. Kamel, and O. Bendjeghaba, “A new Quantum Chaotic Cuckoo search algorithm for data clustering,” Expert Systems with Applications, vol. 96, pp. 358–372, 2018. View at: Publisher Site | Google Scholar
  18. M. Ayoubi and R.-A. Hooshmand, “A new fuzzy optimal allocation of detuned passive filters based on a Nonhomogeneous Cuckoo search algorithm considering resonance constraint,” ISA Transactions, vol. 89, pp. 186–197, 2019. View at: Publisher Site | Google Scholar
  19. M. Z. Ali, N. H. Awad, R. G. Reynolds, and P. N. Suganthan, “A balanced fuzzy Cultural Algorithm with a modified Levy flight search for real parameter optimization,” Information Sciences, vol. 447, pp. 12–35, 2018. View at: Publisher Site | Google Scholar
  20. R. E. Kalman, “A new approach to linear filtering and prediction problems,” Journal of Basic Engineering, vol. 82, no. 1, pp. 35–45, 1960. View at: Publisher Site | Google Scholar
  21. X. Liu, H. Qu, J. Zhao, and P. Yue, “Maximum correntropy square-root cubature Kalman filter with application to SINS/GPS integrated systems,” ISA Transactions, vol. 80, pp. 195–202, 2018. View at: Publisher Site | Google Scholar
  22. S. San Gultekin and J. Paisley, “Nonlinear Kalman filtering with divergence minimization,” IEEE Transactions on Signal Processing, vol. 65, no. 23, pp. 6319–6331, 2017. View at: Publisher Site | Google Scholar
  23. J. E. Solís-Pérez, J. F. Gómez-Aguilar, L. Torres, R. F. Escobar-Jiménez, and J. Reyes-Reyes, “Fitting of experimental data using a fractional Kalman-like observer,” ISA Transactions, vol. 88, pp. 153–169, 2019. View at: Publisher Site | Google Scholar
  24. Y. Sun, X. Wu, J. Cao, Z. Wei, and G. Sun, “Fractional extended Kalman filtering for non-linear fractional system with Lévy noises,” IET Control Theory & Applications, vol. 11, no. 3, pp. 349–358, 2017. View at: Publisher Site | Google Scholar
  25. T. Liu, S. Cheng, Y. Wei, A. Li, and Y. Wang, “Fractional central difference Kalman filter with unknown prior information,” Signal Processing, vol. 154, pp. 294–303, 2019. View at: Publisher Site | Google Scholar
  26. X.-S. Yang, Cuckoo Search and Firefly Algorithm: Theory and Applications, Springer, London, UK, 2014.
  27. Augerat, Known Best Results, [EB/OL], 2013, http://neo.lcc.uma.es/vrp/known-best-results/.
  28. R. Sagban, K. Ruhana Ku-Mahamud, and M. Shahbani Abu Bakar, “Nature-inspired parameter controllers for ACO-based reactive search,” Research Journal of Applied Sciences, Engineering and Technology, vol. 11, no. 1, pp. 109–117, 2015. View at: Publisher Site | Google Scholar
  29. Auto Navi Map, Real-time City Details, [EB/OL], 2019, https://report.amap.com/detail.do?city=350200.

Copyright © 2020 H. Xue. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.


More related articles

 PDF Download Citation Citation
 Download other formatsMore
 Order printed copiesOrder
Views192
Downloads190
Citations

Related articles

Article of the Year Award: Outstanding research contributions of 2020, as selected by our Chief Editors. Read the winning articles.