## 4 Proposed Algorithm
Relative Percentage Deviation (RPD) is used to evaluate the performance of the algorithm. The value of RPD reflects the difference between one solution and the optimal solution in the same instance. The better algorithm corresponds to a smaller PRD value. By observing the change of the value, the difference between the optimization effects of algorithms can be clearly observed. \[\text{RPD}=\frac{Fit(C)-Fit(C)}{Fit(C)}, \tag{13}\]

Equation (13) depicts the calculation of the RPD value. In Equation (13), \(Fit(C\prime)\) represents the optimal fitness for the same instances. ### Parameter Calibration

The initial temperature \(t\) is critical to the performance of the SA algorithm framework. We use a large number of instances to determine the value of \(t\). To the best of our knowledge, no uniform testing benchmark is available for the considered problem. So, we combined some practical scenarios in logistics transportation to generate test instances, such as [22] and [23]. Parameters are designed as follows. (i) We consider the regional parameters with the area radius \(r\in\{1\text{ km},2\text{ km},3\text{ km}\}\) and the number of stations \(m\in[5,10,15]\). (ii) The weight of each package \(w_{i}\) is generated randomly in interval \([0.05,\ g]\), and the parameter \(g\) represents the maximum weight of all packages. The package related parameters are set as \(g\in\{3\text{ kg},4\text{ kg},5\text{ kg}\}\) and \(G\in[80,100,120]\). (iii) The number of UAVs \(h\in\{3,4,5\}\). The special logistic UAV is selected with maximum flight time \(e=0.67\text{ h}\), load capacity \(c=5\text{ kg}\), and average speed \(v=48\text{ km/h}\). There are \(3\times 3\times 3\times 3\times 3=243\) instance combinations. For each instance combination, we compare the algorithm performance of five different values of initial temperature \(t\) of SA framework. Therefore, there are \(3\times 3\times 3\times 3\times 3\times 5=1215\) instances for parameter calibration. One-way analysis of variance technique is used to analyze the experiment results. The mean plots of SA framework with different values of \(t\) are depicted in Figure 3. The smaller RPD values, the better optimization effect. From Figure 3, it can be observed that the differences are statistically significant for \(t\leq 120\), however, they tend to be flat for \(t\geq 120\). So, when \(t\) is set to 120, the algorithm can generate a better solution with fewer iterations. ### Robustness Analysis

In terms of the parameter calibration, \(t=120\) is set to the initial temperature of the SA algorithm framework. To further demonstrate the robustness of the proposed algorithm, we analyze the influence of the five instance parameters on this algorithm. The proposed algorithm is compared with the Local Search (LS) algorithm and the Simulated Annealing (SA) algorithm. As can be seen from Figure 4, the VNS performs better than the other two algorithms. The statistical difference is significant. This is because the VNS algorithm framework can compare the optimal solutions of different neighborhoods, and the final result is closer to the global optimal solution, which can effectively avoid falling into the local optimal situation. Figure 3: The mean plot of the initial temperatures with 95.0 percent Tukey HSD (Honestly Significant Difference) intervals. Interactions between each parameter and the compared algorithms with 95.0 percent Tukey HSD (Honestly Significant Difference) intervals are depicted in Figures 5 and 6. It can be concluded from Figures 5 and 6 that the observed differences are not statistically significant for the proposed algorithm in most cases. Figure 5 shows that the area radius \(r\) has a great effect on the performance of proposed SA algorithm framework. Furthermore, the differences in RPD values are not statistically significant for the number of express cabinet \(m\). This is because the radius of the region directly affects the fitness of the solution.