Locating Movable Parcel Lockers under Stochastic Demands

: The movable unit equipped with a set of lockers has been recently developed as a new mode to improve the e ﬃ ciency of the last mile delivery. Locating a set of movable parcel locker units appropriately is a fundamental factor to promote the merits of movable parcel lockers. However, the di ﬃ culty in determining where to locate movable parcel locker units arises from the stochastic characteristics of demands. Therefore, we propose a robust optimization approach to determine the number of movable parcel locker units and their locations simultaneously with the aim to minimize the operating cost under stochastic demands and mobility restrictions. To reduce the complexity of the optimization model, the non-linear constraints have been transformed into the linear counterparts, resulting in an integer linear programming (ILP) model that can be solved by commercially available mathematical programming solvers. The results from the numerical examples indicate the proposed approach can obtain high robustness with a small extra cost within reasonable time. In addition, it is found that if each unit is equipped with more lockers, fewer movable parcel locker units are required to accommodate the demands with less operating cost, as the demand points can be clustered into a few intensive self-pickup sites.


Introduction
The rapid development of e-commerce has led to great changes in people's shopping mode, and the market share of business-to-consumer (B2C) continues to increase [1,2]. Online shopping transactions of global consumers reached ¥15.57307 trillion in 2015, an increase of 19.9% compared to 2014 [3]. According to Francke and Visser [4], the scale of e-commerce in Europe, North America and Asia Pacific, including the proportion of online consumers, the share of online goods and the amount of online transactions, is considerable. In China, online sales on the day of "double 11" accounted for ¥410.1 billion in 2019, an increase of 30.5% compared to 2018 [5]. The direct impact is that there is a large demand for the last mile delivery.
A large number of orders cannot be efficiently delivered only through home delivery. Courier door-to-door delivery is unable to consolidate fragmented orders [6]. What is more, a research in the U.K. reported that the rate of failed first deliveries of home delivery was 12%, which is caused by the fact that half of U.K. homes are empty from 9 a.m. to 4 p.m. on weekdays largely due to changes in the lifestyle [7,8]. Redelivery is an important source of last-mile delivery costs, and if all packages can be successfully delivered in the first time, it can save ¥7.36933 billion in London [4]. The additional trips have also exacerbated environmental pollution [9,10].
Developing alternative delivery points, called pick-up points (PP) and collection-and-delivery points (CDP), to enable customers to participate in last-mile delivery is a promising solution [11,12]. Iwan et al. [13] classified the alternative delivery points as reception boxes, delivery boxes, controlled access systems, collection points and parcel lockers. The most popular one of these is the parcel locker called the fixed parcel locker in this paper, which supports unattended delivery and solves the problem that the customer is not at home when the courier delivers the goods. Fixed parcel lockers have been launched in many countries around the world [14]. Deutsche Post DHL had more than 13,000 locker service points worldwide by 2019 [15]. InPost in Poland has contributed significantly to the widespread spread of fixed parcel lockers in Europe [16]. The number of fixed lockers in major cities in China has reached 406,000 sets in 2019, an increase of 127,000 sets compared with 2018 [17].
Fixed parcel lockers also have some inherent limitations. First, its location and configuration are fixed, which has poor accessibility and poor ability to respond to change of demands. Moreover, the fixed parcel lockers have high requirements for information and communication technology (ICT) and cannot be set up in old residential areas. In addition, high construction costs are required [18]. For those reasons, a new type of parcel lockers, named movable parcel lockers in this paper, have been recently developed and commercially operated by the providers rather than logistic companies, which means the cost for parcel delivery is not considered by the movable parcel locker providers. For convenience, we refer to a movable parcel locker as a movable unit, on which a number of lockers are amounted. More details on the movable parcel lockers can be obtained from [19].
This paper proposes a robust optimization model with the aim to determine the locations and the number of movable parcel lockers, which are used as the supplement to the fixed parcel lockers, by taking into account of the stochastic characteristics of demands and mobility restrictions, in order to minimize the operating cost for the movable parcel lockers. The rest of this paper is arranged as follows. The next section begins with a description on the problem to be solved and the assumptions made for the model. After that, it focuses on the establishment process of the optimization model proposed with consideration of stochastic demands and the robust optimization method to deal with stochastic demands. In Section 4, a number of numerical experiments were performed to analyze the robustness of solutions and the impacts of key parameters and mobility restrictions on the solutions. Finally, this paper is concluded with a summary of the main research findings and corresponding discussions.

Literature Review
Last-mile delivery methods mainly include door-to-door delivery through vehicles, drones, robots and crowdsourcing strategies, and alternative delivery points where the customers need to go to collect their parcels. Alternative delivery points can be generally divided into supervised pickup points located at local small businesses or convenient shops [16], and unattended dispensing machines or tools (such as parcel lockers and delivery boxes). The door-to-door delivery method has the characteristics of high error rate of delivery addresses, fragmentation and time window mismatch [20], resulting in a high secondary delivery rate and low customer satisfaction. Compared with door-to-door delivery methods, alternative delivery points can facilitate wider service time window, which can significantly alleviate the time window mismatch between customers and couriers. In addition, the unattended delivery point can overcome the limitation of the supervised pickup point being limited by business opening hours, and is a promising solution [16]. In this paper, we will focus on the movable parcel locker which is a new type of the unattended alternative delivery point.
Due to the fact that the type of movable parcel lockers is relatively new, few research works are devoted to it. Therefore, we provide a brief summary on the research works related to the fixed parcel lockers, which can be served as a reference for our work presented in this paper.
Many studies have confirmed that the location of the fixed parcel locker is an important factor influencing its utilization rate and cost. Lemke et al. [21] conducted an online survey of users of InPost's fixed lockers in Poland, and 15% of respondents said that convenient locations could attract them to use lockers more frequently. Kedia et al. [22] found that locker stations with available parking spaces nearby are more attractive to consumers. Lachapelle et al. [23] studied the impact of the locations of fixed lockers on Australian cities. The locations closely integrated with urban design are more appealing.
The location problem has been extensively studied in many fields. Liu et al. [24] designed a bi-level programming model to optimize the location of railway classification yards and used Gurobi software to solve it. Guo et al. [25] studied the integrated optimization of community shuttle stops and routes, and developed a non-dominated sorting genetic based algorithm to solve the proposed bi-objective programing formulation. Guo et al. [26] determined the stop location and running route of the feeder bus service for high-speed rail stations through a bi-objective mathematical model, and used the exact ε-constraint method to solve the problem. However, there is less research on the location optimization of parcel lockers. Deutsch and Golany [20] for the first time quantitatively determined the optimal locations, number and sizes of the fixed lockers with the goal of maximizing total profit. Ji et al. [27] proposed the multi-objective location problem of fixed parcel lockers and aimed to minimize the total cost and energy consumption. Lee et al. [28] established a decision-making system for Korean communities, which can determine the number of locations and the optimal locations of fixed lockers.
Most research on optimizing the locations of parcel lockers assumes that customer demands are deterministic. However, customer demands are stochastic in the real world [29][30][31][32]. Common methods for dealing with stochastic demands include estimating probability distribution [33,34], robust optimization [35,36] and uncertainty theory [37,38]. The traditional classical algorithm or intelligent heuristic algorithm can be used to solve the problem based on the established model and other problem characteristics. Tanonkou et al. [39] established a two-stage non-linear discrete stochastic optimization model to determine the distribution centers, and solved it by a Lagrangian relaxation based approach. Zhang et al. [40] applied a hybrid genetic algorithm to solve the location-routing-inventory problem under stochastic demands. Li et al. [41] designed an improved Benders decomposition algorithm to solve a model with multiple scenarios in the uncertain environment. Wen et al. [31] integrated simplex algorithm, Monte Carlo simulation and genetic algorithm into a hybrid intelligent algorithm in order to solve the capacitated facility location-allocation problem. Brenna et al. [42] proposed a new method to determine the location of electric vehicle charging stations, based on the parallelism between mobility needs and heat equation implemented with finite element method analysis (FEM). This paper assumes that customer demands fluctuate within a certain symmetric interval, and proposes an approach to determine the locations and the number of movable parcel locker units with consideration of stochastic demands. Based on robust optimization, an integer linear programming (ILP) model with the goal of minimizing the operating cost has been established by transforming the non-linear constraints into linear constraints. The model can be solved efficiently using the commercially available mathematical programming solvers (e.g., MATLAB, CPLEX and GUROBI). In this study, we solved the developed optimization models by using the MATLAB solvers. The results obtained from a series of numerical experiments have examined the adaptability and robustness of the proposed approach to the stochastic demands with different level of uncertainty. Finally, a set of experiments with different scales of demand points were performed to analyze the impacts of key parameters on the costs, the locations and the number of movable parcel lockers.

Robust Optimization Model under Stochastic Demands
This section begins with a description on the problem and the assumptions needed to develop the optimization model. Before introducing the robust optimization model for stochastic demands, the deterministic version is presented with the aim to minimize the operating cost. Finally, the robust optimization model is transformed into an equivalent ILP model. In the integer programming (IP) model, some or all of the decision variables must take integer values. If the decision variables can only take 0 or 1, the model is called a 0-1 IP model. Except for integer limitation, if the objective function and all constraints are linear, the model is an ILP model. The cutting plane method and the branch and bound method are frequently used to solve the IP problem. The cutting plane method can find the optimal solution of the IP problem at an integer pole of the reduced feasible region by adding a finite number of cutting plane constraints. The branch and bound method finds the optimal solution by repeatedly dividing the solution space into smaller and smaller subsets, and calculating the upper or lower bound of the solutions in each subset [43,44]. The decision variables of the ILP model developed in this paper include integer and 0-1 variables, and the model is solved by the intlinprog (MATLAB) [45].

Problem Description and Assumptions
At present, fixed parcel lockers are already deployed in many sites (e.g., residential communities and shopping centers). However, these fixed lockers are less flexible to response to the demand variation. For instance, it has been experienced that the fixed parcel lockers are not adequate to meet the increased demands during novel coronavirus pandemia. Therefore, this paper attempts to propose an approach to determining the number of movable parcel locker units and the location for each unit as supplement to the fixed lockers with consideration of the stochastic characteristic of demands.
Without loss of generality, this paper assumes a network G(I, E) consisting of I, the set of demand points and E, the set of edges connecting two demand points (incomplete connection). Note that the demand point in this paper refers to a community or a part of a community whose demand can only be allocated to one self-pickup site which is the nearest one. The communities discussed in this paper are those where the union of the time windows required by all customers occupies almost the entire daytime, so that there is no or small service gap for the movable parcel locker units in these communities during the daytime. Therefore, it is difficult to dynamically operate the movable parcel locker units in such circumstances, as it takes the entire daytime to empty all parcels in a unit, and also frequent movement may result in a reduction in the profit, due to potential increase in the relocating cost from long-term point of view. Furthermore, a survey of the movable parcel locker unit operator implies that daily-based operation is a favorable mode as it can alleviate extensive relocation and reduce unnecessary congestions or accidents and labor force involved in the relocating task during daytime. On the other hand, some operators worry about the unexpected movement or damage in unattended movable parcel locker units parked at some points rather than the depot. For these reasons, we assume the movable parcel locker units are operated in such a way that they are dispatched from the depot in the morning and returned back in the evening for easy replenishment and maintenance. However, such a periodic operation mode should solve the issue that some parcels may still remain until the time when movable parcel locker units should be returned back. There are two possible ways to deal with this issue. If fixed parcel lockers nearby are available, the remaining parcels will be transferred to the fixed parcel lockers. Otherwise, the remaining parcels will be taken back to the depot and redistributed as part of demands in the following day. According to the above discussion, this paper summarizes this operation problem of the movable parcel locker units as follows. A number of movable parcel locker units are dispatched from the depot to a set of demand sites called self-pickup sites in this paper (i.e., a subset of whole demand points) every morning, and are returned back to the depot for maintenance and power replenishment at night. Therefore, it is necessary to decide which demand points should be chosen to accommodate the movable parcel locker units according to the daily demand with minimum cost. In detail, it is necessary to find out some areas called self-pickup sites from the demand point set to place movable parcel locker units, and determine the number of units in each area. It is assumed that customers are most willing to pick up or return the goods from or to the nearest self-pickup site. We aim to minimize the operating cost including the rent for land, investment cost, maintenance cost and travel cost with the limitation of stochastic demands.
An example is given to illustrate our problem. In Figure 1, "SPKS" and "DP" represent self-pickup sites and demand points, respectively. As shown in Figure 1, there are 25 demand points in a distribution area, and the demand points may be different every day, which is distinguished by different colors. In Beijing, the scale of the distribution area is about 3-6 km 2 . Figure 1 shows the locations and service objects of movable parcel locker units for three days. The circles represent demand points, and the filled circles represent self-pickup sites that are selected from those demand points and are locations of movable parcel locker units. In Figure 1, each dashed line indicates which self-pickup site each demand point is attached to. Each movable parcel locker unit stays at one self-pickup site to meet all demands within its radiation range. The radiation range depends on the maximum walking distance acceptable to the customers. For example, demand point 3 is one of the self-pickup sites on the first day, which covers demand points 1, 2 and 3. Therefore, the customers of demand points 1, 2 and 3 should pick up or return parcels from or to the movable parcel locker unit at demand point 3. This example assumes that one movable parcel locker unit can satisfy requests from all nearby demand points. Otherwise, multiple units need to be placed. Figure 1 shows that three movable parcel locker units meet the three-day demands of all demand points in the distribution area. We can find that nine fixed parcel lockers need to be set up under the same demands based on the fact that once the fixed lockers are set up, it is difficult to change their locations.

Optimization Model under Stochastic Demands
Before presenting the optimal model, the symbols involved in the model are shown in Table 1.

Optimization Model under Deterministic Demands
The objective function with the aim to minimize the operating cost can be expressed as where x i and A represent the total number of lockers set up at the self-pickup site i, in order to meet the delivery demands of itself and all the demand points allocated to it as much as possible under stochastic demands, and the number of lockers equipped on a movable parcel locker unit, respectively. Thus, "n i = x i /A ≥ 0" indicates the number of movable parcel locker units at the self-pickup site i. In fact, only when n i or x i > 0, the demand point i is chosen as a self-pickup site to accommodate movable parcel locker units. Note that the operation to round up to nearest integer for the number of movable parcel locker units (i.e., x i /A ) is not included in Equation (1) in order to simplify the calculation. As shown in Equation (1), the operating cost of movable parcel locker units consists of two parts. The first part, i∈I c f x i /A, represents the purchase and maintenance costs of movable parcel locker units. The second part, i∈I c s i x i /A, shows the cost for the land use and transportation of movable parcel locker units, with c s i denoting the land rental rate for a movable parcel locker unit at site i and the transportation cost of a unit from the depot to site i. If some road sections are closed to vehicles at nights and this causes the movable parcel locker units to detour back to the depot from some sites, the corresponding sites have a larger c s i , which subsequently reduces the probability of those sites being selected as self-pickup sites and increases the cost compared to the normal situation. Clearly, the cost of the first part depends on the number of movable parcel locker units, while the cost of the second part is also related to the location of the self-pickup sites selected.
Then, the two unit costs are combined as a generalized cost of each locker, which is shown in Equation (2).
Equations (1) and (2) attempt to minimize the operating cost, including the land rent, investment, maintenance and transportation cost, with respect to the number of lockers and their locations (i.e., x i ).
A set of constraints have been formulated and shown as (3)-(10).
Constraints (3)(4)(5) show the value range of x i at each site. If x i = 0, then no movable parcel locker unit is placed at the demand point i. When x i > 0, it means that the demand point i is selected as a self-pickup site to accommodate a number of movable parcel locker units. Therefore (6), y ij is 1, if the jth demand point is allocated to the ith self-pickup site, which means the customers living in the jth demand point have to go to the ith self-pickup site to collect or sent parcels, and 0, otherwise. Constraint (7) ensures that a demand point can only be assigned to one self-pickup site. Constraints (8) and (9) mean that if movable locker units are placed at the ith demand point, the demand of the ith demand point must be served by those units. Note that W in Constraint (8) is a big positive number. Constraint (10) indicates that if the jth demand point is assigned to the ith site, then some lockers must be set up at the ith demand point. In Constraint (11), l ij represents the distance between the ith self-pickup site and the jth demand point. Since it is assumed that customers are most willing to pick up or return the goods from or to the nearest self-pickup site, Constraint (11) means that each demand point is allocated to one self-pickup site which is the nearest one. Constraint (12) restricts the distance from the jth demand point to its own self-pickup site i not to exceed the maximum walking distance r acceptable to customers. In Constraint (13), d j represents the demand at the jth demand point except for the demand satisfied by fixed lockers. Constraint (13) indicates that the number of lockers provided at the ith self-pickup site needs to meet all the demands in the surrounding sites that have been allocated to the ith site.
Obviously, the Constraint (11) is a non-linear constraint due to the min function. To simplify the calculation process, the Constraint (11) is converted to the linear Counterparts (14) and (15). My where M is a very large positive number and should be large enough to satisfy "max i∈I/I indicates that demand point i is the kth closest to demand point j. y j(1) ij = 1 means that the demand point j is allocated to the demand point i closest to it. That is, the demand point closest to it is selected as a self-pickup site. Correspondingly, the number of lockers set up at the self-pickup site i is x j(1) i . Otherwise, the demand point j is allocated to the demand point that is the kth closest to it.
Equations (14) and (15) indicate that each demand point is allocated to the self-pickup site closest to it. The proof process of Equations (14) and (15) is as follows.

Lemma 2.
From i∈I y ij = 1 and y j(k) ij = 1 , it can be seen that y j(s) Lemma 2 shows that the jth demand point must be allocated to only one self-pickup site.
Proof. (14) and (15) will be proved by the following scenarios, based on the Lemma 1 and 2. • k > 1 According to Lemma 2, y j(k) ij and y j(k) ij = 1. In the second scenario, the demand point i is selected as a self-pickup site, which is the kth closest site to the jth demand point, and the jth demand point is assigned to the demand point i with the kth shortest distance to it.
In the third scenario, the demand point i is selected as a self-pickup site, which is the kth closest site to the jth demand point. However, the jth demand point is assigned to the demand point with the hth shortest distance to it.
Whether k > h or k < h, since y j(h) ij = 1 and i∈I y ij = 1, y j(k) It can be seen from the above proof that the non-linear Constraint (11) is equivalent to its linear Counterparts (14) and (15). This means the solution space is unlikely affected by such transformation and therefore the optimal solution can be guaranteed. Furthermore, a literature review on the works related to the transformation of the non-linear constraints into the linear counterparts does not show any implication of such transformation on the model precision and robustness. On the other hand, it has been frequently reported that the easiness to derive the optimal solution after the transformation of the non-linear constraints.

Robust Optimization
It is assumed that the demand that is not satisfied by the fixed lockers at each demand point is a random variable. However, it is assumed that the distribution of the random variable is unknown, but we have knowledge of the range within which it varies. For example, the unmet demand d j of the jth demand point varies in the symmetric interval d  − ∆d j , d  + ∆d j . According to the robust optimization method proposed by Bertsimas and Sim [35], Constraint (13) in the model under deterministic demand can be expressed as where J i represents the set of demand points assigned to the ith self-pickup site, S i ⊆ J i the set of demand points in J i whose demand variation is ∆d j , Γ i a parameter, which can take any value in the range [0, |J i |] (not necessarily integer) and t i the demand point in J i whose demand variation is less than the largest one ∆d j . Note that we can use Γ i to adjust the trade-off between robustness and risk. When Γ i = |J i |, the demand of all demand points is the largest, which represents the most conservative situation [35]. However, it has little chance that demands of all the demand points take the maximum variation at the same time.
For certain y * i and For certain y * i and Γ i , Equation (17) attempts to maximize the total demand variation at self-pickup site i with respect to the set J i of the demand points allocated to self-pickup site i.
Dual of problem of (17)- (19) can be expressed as follows: s.t. p ij + q i ≥ ∆d j y * ij , ∀i ∈ I/I , j ∈ I p ij ≥ 0, ∀i ∈ I/I , j ∈ I q i ≥ 0, ∀i ∈ I/I So far, the stochastic constraint has been transformed into deterministic constraints. Therefore, the equivalent deterministic model of the stochastic model can be summarized as an ILP model as follows. min i∈I c i x i x i ≥ y ij , ∀ i, j ∈ I/I , i = j (31) x i ≥ y ij , ∀ i ∈ I/I , j ∈ I (32) l ij y ij ≤ r, ∀ i ∈ I/I , j ∈ I p ij + q i ≥ ∆d j y ij , ∀i ∈ I/I , j ∈ I p ij ≥ 0, ∀i ∈ I/I , j ∈ I The above ILP model can be solved using a number of commercially available mathematical programming solvers. In this paper, we adopted the solver, intlinprog (MATLAB) [45], to solve the developed model. Figure 2 shows the structure for our robust optimization algorithm to determine the locations and the number of movable parcel lockers simultaneously. First, according to the connection status and mobility restriction information of the road sections in the distribution network, the shortest path tree from the depot to all demand points is obtained by using the shortest path algorithms (e.g., Dijkstra algorithm used by this work). Then, based on the shortest path tree obtained, the shortest distances between each pair of demand points and from the depot to each demand point can be derived.
Subsequently, c s i and c i can be estimated based on the shortest distance from the depot to each demand point. Finally, all parameters are input into the stochastic optimization model to obtain robust solutions to the optimization problem.

Experiments and Results
To evaluate the proposed approach to locating movable parcel locker units under stochastic demands, a series of numerical experiments have been performed. This section begins with an introduction to the parameter settings for the experiments presented in this paper. The robust results with different values of Γ i and probability bounds of constraint violation are presented and discussed in Section 4.2. A set of experiments with different scales of demand points are reported, in Section 4.3, with a discussion on the impacts of the maximum walking distance and the sizes of the movable parcel locker units on the optimization results. Finally, the impacts of mobility restrictions on the optimization results are analyzed in Section 4.4. All experiments were implemented in the MATLAB environment [45] on an x64-PC with an Intel Core i7-8550U 1.80 GHz CPU and 8.0 GB of RAM.

Parameter Settings
First, the costs of a movable parcel locker unit are estimated based on the data (as shown in Table 2) obtained from a survey to a movable parcel locker unit provider. At present, there are two types of movable parcel locker units emerged in the Chinese market. One is the movable parcel locker unit equipped with 64 lockers (64 MPLs), and the other has 128 lockers (128 MPLs). Their carriers are electric vehicles equipped normally with a battery pack consisting of two battery of 1.2 kWh. When a movable parcel locker unit is parked, the average power of the system does not exceed 40 W, and the average power during driving does not exceed 400 W [19]. It is assumed that the movable parcel locker unit works 12 h a day, and the average travel time is about 2 h. In addition, each unit is equipped with a solar panel to supplement electrical energy. Therefore, it is conservatively believed that the movable parcel locker unit can be used for 1 day when fully charged. The charging process of the movable parcel locker units can be completed in the residential areas [46] or the depot. This paper assumes that the movable parcel locker units return to the depot every night for charging and maintenance. The costs of the two types of movable parcel locker units are detailed below. Taking the Chinese market as an example, the price of 64 MPLs and 128 MPLs is ¥20,000 and ¥32,000, respectively. Therefore, the purchase cost is ¥5.48/day (64 MPLs) or ¥8.77/day (128 MPLs) assuming a service life of 10 years.

Maintenance Cost
The maintenance cost includes battery replacement cost and communication cost. A battery pack containing two batteries is replaced every two years. The price of each battery is ¥400, so the battery replacement cost is ¥1.1/day. With the assistance of cell phone APPs (e.g., WeChat), collecting or sending parcels can be accomplished in unattended manner, but there is a cost for the communication.
Although it seems that the communication cost should depend on the number of parcels to deliver, the movable parcel locker units can often obtain special offer from mobile network operators. Therefore, we assume the communication costs of 64 MPLs and 128 MPLs are ¥2.56/day and ¥5.12/day on average, respectively. Therefore, the maintenance cost is ¥3.66/day (64 MPLs) or ¥6.22/day (128 MPLs).

Travel Cost
The travel cost includes charging cost and labor cost, which may be different for different demand points due to possible different travel time to them. According to the previous assumptions, a movable parcel locker unit is moved by the driver for two hours on average every day and needs to be charged once a day. Therefore, two batteries need to consume 3 kWh of electricity (charging efficiency is 80%). The flat peak charge standard of general industry and commerce in the urban area of Beijing is ¥0.78/kWh [47], so the charging cost is ¥1.17/hour. Taking Beijing as an example, the average hourly salary of staff in private enterprises is about ¥17.56 [48]. Therefore, the travel cost is ¥18.73/hour. Traditional courier vehicles need to pay about ¥300/month for special cases (mainly referring to illegal parking). Therefore, in this paper, the rent of a movable parcel locker unit is set as about ¥300/month, or considering the sizes of the units, set as about ¥3/m 2 /day. For different demand points, the rent can fluctuate around this value.
The stochastic demands of customers are generated randomly with the average demand d  and the fluctuating demand ∆d j . The set of maximum walking distance was (150 m, 300 m and 500 m). The number of demand points was assumed to be one of (10, 30, 50, 100 and 150). The value of Γ i is in the range [0,|J i |]. In each experiment, in order to analyze the impacts of one parameter on the results, the other parameters were kept the same. It was clear that the demands were deterministic, when Γ i = 0. On the other hand, "Γ i = |J i |" means that the demand variation of all demand points reaches the maximum value.

The Robustness of Solutions
In this section, we set the number of demand points to 50, the maximum walking distance to 150 m, the value of average demand d j being in the range (from 50 to 150) and the value of maximum demand variation ∆d j being in the range (from 10 to 40). Then, 64 MPLs were used to complete delivery tasks. Robustness is an important criterion for evaluating the performance of an optimization technique [49]. Therefore, we performed a series of numerical experiments to analyze the robustness of solutions.
Γ is the key parameter to control the robustness for the optimization problem of movable parcel locker units under stochastic demands. The larger value of Γ i means the higher level of protection, but the larger gap between its corresponding objective function value and the optimal objective value for deterministic demand. Figure 3 reflects the impact of Γ i on the optimal cost calculated by Equation (2). "Relative cost ratio" indicates the relative gap between the optimal cost with different values of Γ i and the optimal cost under deterministic demand. It can be found that the optimal cost does not increase significantly with the increase of Γ i (0 < Γ i < |J i |). When Γ i = 5, the optimal cost increases by 20.76%. When Γ i = 45, the optimal cost increases by 24.70%. Figure 4 shows the optimal cost with different probability bounds of constraint violation given in the equation of Bound 3 in reference [35] and reflecting the robust levels. As shown in Figure 4, the probability bound of constraint violation significantly decreases with the small increase in cost, reflecting robustness greatly increase.   Table 3 presents the results of a set of experiments, showing the trade-off between optimal cost and robustness with different values of Γ i and the computational performance of the proposed approach. The column "Lockers" indicates the optimal number of lockers provided in each numerical experiment with consideration of stochastic demands. The column "Relative locker ratio" means the relative gap between the optimal number of lockers under the stochastic demands with different level of uncertainty and the optimal number of lockers under deterministic demand. The column "Optgap" means the relative gap to the optimal solutions. The column "Time" represents the total running time to find optimal solutions. As shown in Table 3, when Γ i = 18, i.e., only 36% of the number of uncertain data, the probability bound of constraint violation is less than 1%. Moreover, the larger values of Γ i can obtain higher robustness in small increments of relative cost ratio and relative locker ratio. The optimal gap of 0 means that the reliability of optimal solutions is too high to display due to the limit of the display resolution. All optimal solutions were obtained within a minute. Based on the experimental results, the following conclusions can be drawn. The robust optimization method developed in this paper is less conservative. That is to say, the high robustness can be obtained without excessively increasing the value of Γ i . On the other hand, the transformation of the original non-linear problem into the linear counterpart reduces the problem complexity, so that no specific algorithm is required to develop, and the solution generated by the existing solver can be obtained within a reasonable time.

The Impacts of Key Parameters on the Optimization Results
This paper is devoted to locating the movable parcel locker units under stochastic demands and proposes an approach to determine the locations and the number of movable parcel locker units with the goal of minimizing the operating cost. In Section 4.2, the impacts of Γ i on the robustness of the solutions are analyzed. This section attempts to discuss changes in optimization results when the values of other parameters change. In detail, this section presents a series of experiments with five scales of demand points (i.e., (10,30,50, 100 and 150)) designed to analyze the impacts of the size of movable parcel locker unit and the maximum walking distance on the costs, the locations and the number of movable parcel locker units. For this purpose, the following parameters are employed. The sizes of movable parcel locker units examined in this work are 64 and 128 MPLs, while three different maximum walking distances (i.e., (150 m, 300 m and 500 m)) are used for the comparison purpose. Note that Γ i is set to 0.7*|J i | for all experiments presented in this section. Figure 5 shows the experimental results. The abscissa represents the combination of the number of demand points, the size of movable parcel locker unit and the maximum walking distance. For example, (50, 64 and 150) represents a scenario in which the number of demand points is 50, the size of a movable parcel locker unit 64 and the maximum walking distance 150 m. Figure 5a shows the relative gap in the number of self-pickup sites for different maximum walking distances under the same demand point number and the size of movable parcel locker unit. Figure 5b- Figure 5a shows the impacts of each parameter on the number of self-pickup sites. The impacts of one parameter on the results were analyzed by keeping the other two parameters constant. As can be seen from Figure 5a, only the maximum walking distance had a significant impact on the number of self-pickup sites. With the increase of the maximum walking distance, the number of self-pickup sites decreased considerably. Figure 5b shows the impacts of each parameter on the number of movable parcel locker units. As shown in Figure 5b, the same fleet size of 64 MPLs and 128 MPLs was required by less demand points. However, the fleet size of 128 MPLs was smaller than that of 64 MPLs for the demands of more demand points. Similarly, it was found that increasing the maximum walking distance could reduce the number of movable parcel locker units. These imply that if each unit was equipped with more lockers, fewer movable parcel locker units were required to accommodate the demands, as the demand points can be clustered into a few intensive self-pickup sites.

Impacts on the Costs
The impacts of parameters on the generalized cost calculated by Equation (2) and the real cost calculated by Equation (1) based on n i are shown in Figure 5c,d. The generalized and real costs apply to cases with a locker as a unit and a movable parcel locker as a unit, respectively. As shown in Figure 5c, the generalized cost required by 128 MPLs was much less than that required by 64 MPLs, and the relative gap was about 44%. This is because the total generalized cost depends on the generalized cost of a locker under the same demands. In this paper, a locker of 128 MPLs saved 44% of the generalized cost compared to 64 MPLs (for details, see the discussion in Section 4.1). The impacts of parameters on the real cost were similar to their impacts on the number of movable parcel locker units. Furthermore, 64 MPLs can respond to the demands of a small number of demand points with less real cost. However, 128 MPLs can save the real cost at a larger demand point scale. This was because the real cost depends heavily on the number of movable parcel locker units. As the number of demand points and the maximum walking distance increase, the quantity advantage of 128 MPLs makes up for its disadvantage in unit cost. Therefore, the movable parcel locker units equipped with more lockers can accommodate the demands with less real cost, as the demand points can be clustered into a few intensive self-pickup sites. However, it can always save the generalized cost to dispatch the idle lockers from 128 MPLs.

The Impacts of Mobility Restrictions on the Optimization Results
Due to the road maintenance and other activities at nights, some road sections are closed to vehicles at nights (i.e., section mobility restriction). In addition, a number of demand points cannot be selected as self-pickup sites due to mobility restrictions such as forbidding the parking of the movable parcel locker units. On the other hand, some demand points must be selected as self-pickup sites due to mobility restrictions such as their customers' inability or inconvenience to go to other sites to pick up or return parcels (e.g., the aged residential communities). In this work, we regard those two situations as "site mobility restriction".
Based on the above situations, this section briefly analyzes the impacts of mobility restrictions on the optimization results. In Figure 6, {12} represents the depot, and the remaining nodes are demand points. The demand value range of each demand point and the Γ i value are the same as those in Section 4.3. The maximum walking distance was 500 m. Figure 6a shows the link connectivity in the distribution network under normal situation, while the shortest paths between demand points under normal situation are depicted in Figure 6b. In Figure 6c, vehicles are not allowed to drive on arcs (2, 6), (3,12) and (4,11) at night, demand points 4 and 9 forbid the parking of express vehicles, and customers at demand point 1 cannot pick up or return goods from other sites. Under the section mobility restriction illustrated in Figure 6c, the shortest paths between demand points at night are shown in Figure 6d. Note that since it is assumed that there are no road sections restricting vehicle movement during the day, the shortest paths between demand points during the day are the same as those under the normal situation. As shown in Table 4, based on the approach proposed in this paper, the self-pickup sites, the demand points covered by each site and the number of lockers at each site are determined before starting the work of the day. For example, under the normal situation, demand point 2 is selected as the self-pickup site covering demand points 2 and 6, and 244 lockers are set up at the self-pickup site 2. Under the site mobility restriction, demand point 1 is selected as the self-pickup site, and demand points 4 and 9 are not selected as the self-pickup sites, which is different from that under the normal situation. Therefore, the site mobility restriction in this example will increase the cost. Otherwise, if the site mobility restriction does not affect the self-pickup site scheme, the site mobility restriction will not increase the cost. In addition, because it is forbidden to traverse on some arcs at night and this causes the movable parcel locker units to detour back to the depot from some self-pickup sites, as compared with the normal situation. As a result, the section mobility restriction in this example will increase the cost. In this example, the cost under mobility restrictions has increased by 6.85%, as compared to the normal situation.

Conclusions
The work presented in this paper aims to develop an optimization model to locate movable parcel locker units and determine the number of units under stochastic demands. An ILP model has been established based on a robust optimization approach. Γ plays a key role in coordinating robustness and risk. A series of numerical experiments have been performed to examine the adaptability and robustness of the solutions with different values of Γ i and prove that the robust optimization approach adopted in this paper is less conservative. Furthermore, the original model with the non-linear characteristics have been converted into an ILP model, which can be solved by the existing solvers, and no additional effort is required to construct a specific algorithm to solve the model.
The size of movable parcel locker unit and the maximum walking distance are the key parameters affecting the optimization results. The impacts of these parameters on the number of movable parcel locker units, the total generalized cost, the total real cost and the number of self-pickup sites are analyzed based on the results obtained from a set of experiments with different demand point scales. The following conclusions can be drawn. The movable parcel locker units with fewer lockers are suitable for meeting a small number of scattered demands. However, if each unit is equipped with more lockers, fewer movable parcel locker units are required to accommodate the demands with less operating cost, as the demand points can be clustered into a few intensive self-pickup sites. In addition, the number of self-pickup sites can be significantly reduced by increasing the maximum walking distance.
Finally, this paper analyzes the impacts of two types of mobility restrictions, namely site mobility restrictions and section mobility restrictions, on the optimization results. The results show that only when the site mobility restriction affects the self-pickup site scheme of the normal situation, the operating cost will be increased. On the other hand, under the section mobility restriction, if the movable parcel locker units have to take a detour to return to the depot from some self-pickup sites, compared with the normal situation. Then, the section mobility restriction will increase the cost.
As a new concept of last mile delivery, movable parcel lockers have not yet gained widespread popularity. The research in this paper attempts to provide a theoretical support for decision-makers to use movable parcel lockers and formulate implementation plans.