Robust depot location and inventory management in bike-sharing systems under government supervision

Abstract

Bike-sharing systems offer an efficient and accessible option for short-distance transportation. However, the boom of bike-sharing systems also introduced challenges to urban transportation systems, such as congestion due to oversupplying bikes and low service quality caused by malfunctioning bikes. To address these challenges, government supervision has been conducted to improve operational orders in bike-sharing systems. For instance, the government may restrict the bike inventory in certain regions to avoid congestion caused by oversupplying bikes, and supervise the operational condition of the bike fleet to improve the satisfaction of users. Therefore, bike-sharing system operators should incorporate the influence of government supervision when formulating related decisions. This study introduces an optimization model for bike-sharing systems that simultaneously addresses depot location, inventory management, and bike rebalancing decisions under government supervision. Specifically, we consider government supervision aimed at controlling bike oversupply and managing the maintenance of malfunctioning bikes. Moreover, to address the uncertainty of travel demand and government supervision strategy, a distributionally robust optimization (DRO) model that incorporates both moment-based and sample-based information is developed and reformulated. The numerical results imply that the DRO model surpasses benchmark models. Furthermore, effect analyses are conducted, and we provide valuable managerial insights for government supervisors and bike-sharing system operators. For government supervisors, conducting the supervision could effectively mitigate the oversupply of bikes. While for bike-sharing system operators, we notice that the bike allocation restriction and malfunctioning bikes may significantly affect their profits. Accordingly, they should optimize the bike rebalancing operations to avoid such penalties and improve the operational profits.

Graphical abstract

Keywords

bike-sharing system / supervision strategy / depot location / inventory management / distributionally robust optimization

1 Introduction

Bike-sharing systems offer travelers a convenient option for short-distance transportation. Due to their convenience and high accessibility, bike-sharing systems have attracted millions of users globally. Their rapid expansion has also played a role in reducing carbon emissions (Qin et al., 2023). According to the annual report of Hellobike, one of China’s largest free-floating bike-sharing systems, the platform helped cut carbon emissions by over 1.9 million tons in 2023.

However, the boom in the bike-sharing market may also cause troubles in urban transportation systems, particularly for free-floating bike-sharing services. Since there are no restrictions on bike stations and racks in such systems, users could park bikes arbitrarily at any permitted locations, and there are no restrictions on the number of parked bikes. If the bike supply in a certain region is too high, the oversupply of bikes may occupy the roads and cause traffic congestion. For instance, on April 4, 2017, approximately 10,000 bikes flooded Shenzhen Bay Park in Shenzhen, China, leading to severe congestion (He et al., 2020). Moreover, Fig. 1(a) illustrates how a large number of parked bikes on the roadside disrupted traffic. Another issue for the operator is bike malfunctions (Kaspi et al., 2016), because some core components may deteriorate over time, affecting the availability of bikes. For instance, Fig. 1(b) depicts malfunctioning bikes abandoned on the roadside. These unusable bikes reduce service quality, as users may become dissatisfied upon encountering malfunctioning bikes (Xiao and Wang, 2020). These operational troubles may bring challenges to urban transportation systems and hinder the continued growth of bike-sharing services.

To address these operational challenges in bike-sharing systems, many cities and governments have implemented supervision for bike-sharing systems, which contributes to improving market operations and regulatory compliance. For instance, in Singapore, the government restricts the supply of bikes to mitigate congestion caused by oversupply. New bike-sharing operators are required to obtain a “sandbox license”, which limits the maximum of bikes they can operate. Moreover, operators must allocate sufficient resources to manage their fleets and address issues of illegally parked bikes. Since the introduction of these supervision strategies in 2019, the incidence of illegal parking has decreased significantly (Ministry of Transport in Singapore, 2022). Similarly, in China, many cities have implemented regulatory measures to control bike oversupply and illegal parking. Operators with poor operational performance face penalties and may even be banned from further operation.

Though government supervision helps enhance efficiency and social welfare, conducting supervision in all operational regions in large cities during the whole day may require a significant supervision cost. Instead, the supervisor can choose only a proportion of the operational regions to conduct supervision during some operational periods (Zhao et al., 2022). In large cities in China, such as Shenzhen and Xiamen, the transportation managers have implemented such a “random supervision strategy", i.e., to choose the supervised region randomly, to effectively manage the operational order of bike-sharing systems (Shenzhen Municipal Transportation Bureau, 2017; Xiamen Municipal Law Enforcement Bureau of Urban Management, 2025). By conducting such random supervision strategies, the supervisor can effectively address the congestion caused by oversupplying bikes with less cost spending on the supervision. In addition to the supervision of oversupplying bikes, the supervision of malfunctioning bikes is also frequently conducted to address the concern of malfunctioning bikes and improve the service quality. Accordingly, we consider the supervision for both oversupplying and malfunctioning bikes in this paper.

From the perspective of the operators, if the supervisor detects oversupplying and malfunctioning bikes, they are usually punished with extra costs. Accordingly, rather than following traditional operational strategies that aim to optimize more users’ demand and minimize operational costs, the operators should reconsider their operational strategies under the supervision of the government to avoid extra penalty costs. For example, the allocation of bikes among different regions will be adjusted and restricted to avoid oversupplying and the corresponding penalty cost. To address the concerns of oversupplying and malfunctioning bikes, many bike-sharing systems have their own depots located in the operational regions. The depots in bike-sharing systems can be used to accommodate malfunctioning and oversupplying bikes to avoid punishment by the supervisor. In addition, bike rebalancing operations are frequently implemented to address the imbalanced distribution of travel demands and bikes and move the malfunctioning and oversupplying bikes to depots or other regions to avoid punishment.

Moreover, various uncertainties also challenge the operation of bike-sharing systems. One key concern for operators is demand uncertainty, which can impact bike allocation and rebalancing decisions. However, operators may never fully acquire the true distribution of travel demand. In addition to the uncertainty of travel demand, the government’s supervision strategy is unknown to the operators due to the randomly selected supervised regions. Accordingly, the operators should incorporate both the uncertainty of travel demand and supervision strategy when making related decisions.

To address these concerns, this study focuses on strategic and operational challenges within bike-sharing systems, considering the supervision strategy conducted by the government supervisor. Two kinds of supervision, i.e., supervision of malfunctioning bikes and oversupplying bikes, are incorporated into this study. Considering the uncertainties in demand and government supervision strategies, we develop a stochastic optimization model. Furthermore, we extend this model as a two-stage distributionally robust optimization (DRO) model. Before the realization of travel demand and supervision strategies, we determine the strategic decisions like depot locations and allocate bikes across different regions in the first stage. Then, after demand and supervision strategies are realized, we optimize the operational decisions like bike rebalancing and malfunctioning bike transportation between depots and regions throughout daytime operations in the second stage. The proposed DRO model can be further reformulated into a tractable form. We summarize the key contributions of our study as follows:

First, as far as we know, this study is the first to optimize the strategic and operational decisions considering the supervision strategy in bike-sharing systems. Though several papers (Zhao et al., 2022) considered the government supervision and its effect on the bike-sharing systems, they focused on the effect analysis of the government supervision rather than the operation strategy of the operator under supervision. In addition, the supervision strategy uncertainty is also ignored in existing research. Thus, we jointly optimize the strategic and operational decisions in this paper, considering the uncertainty of the supervision strategy conducted by government supervisors. We aim to provide managerial insights for both the supervisors and the operators of bike-sharing systems to improve the operational order of bike-sharing systems.

Second, various kinds of uncertainties, such as the uncertain travel demand and the supervision strategy, are incorporated into the model. To tackle the challenges of uncertainty and data limitation, we develop a DRO model that effectively benefits from the limited available data. We construct an ambiguity set incorporating both moment information and sample-based information from historical data, then reformulate it as a computationally tractable form.

Third, numerical experiments using real-world operational data imply that the DRO model surpasses several benchmark models in out-of-sample testing. In addition, we conduct effect analyses and derive managerial insights on operational strategies for operators of bike-sharing systems considering government supervision.

The structure of this paper is as follows: Section 2 reviews relevant literature, including strategic and operational challenges in bike-sharing systems, government supervision, and optimization methods for handling uncertainty. Section 3 defines the research problem and formulates an optimization model that incorporates supervision in bike-sharing systems. Section 4 addresses demand and supervision uncertainties by introducing a stochastic optimization model, which is further extended into a two-stage DRO framework using a novel ambiguity set. Section 5 presents numerical experiment results and provides managerial insights for both supervisors and bike-sharing system operators. Finally, Section 6 summarizes the study, with proofs of propositions included in the Appendix.

2 Literature review

In this section, we provide a comprehensive review of the literature related to this study. Specifically, we focus on three aspects, including the strategic and operational problems in bike-sharing systems, the supervision in bike-sharing systems, and the optimization approaches for uncertainty in Section 2.1 to 2.3. In addition, we summarize the previous literature and conclude the research gaps in Section 2.4.

2.1 Strategic and operational problems in bike-sharing systems

Decisions in bike-sharing systems span strategic, tactical, and operational levels. At the strategic level, key decisions involve determining bike station locations and designing the overall network. In addition to the strategic level decisions, tactical and operational decisions encompass bike fleet planning (Shu et al., 2013), inventory management (Schuijbroek et al., 2017), and bike rebalancing operations (Raviv et al., 2013). For a comprehensive review of related research in bike-sharing systems, interested readers can refer to Fishman (2016), and Shui and Szeto (2020).

Strategic decisions focus on optimizing facility locations, including station location and bike lane planning in docked bike-sharing systems, to enhance user convenience. Lin and Yang (2011) firstly discussed the station location problem to optimize the total system cost. Furthermore, Lin et al. (2013) integrated inventory management and bike-sharing network design, jointly optimizing strategic and tactical decisions. Frade and Ribeiro (2015) developed a framework to optimize station locations, aiming to maximize service levels within a given budget. Additionally, they considered both station location and capacity design, employing a queuing model to evaluate service levels. Liu et al. (2022) introduced a bike lane planning optimization approach leveraging bike trajectory data.

In addition, the strategic decisions are also jointly optimized with the operational decisions. Nikiforiadis et al. (2021) introduced a joint optimization model for station location and rebalancing operation. They introduced a multi-objective optimization model to satisfy more demand while minimizing the requirement for bike rebalancing. Similarly, Fu et al. (2022b) optimized station location and rebalancing vehicle deployment using a robust optimization model that accounts for uncertainty. They proposed a constraint-generation approach to solve the proposed model. Song et al. (2024) introduced a bi-level station location model for both bikes and e-bikes, incorporating users’ route choices and roaming delays, and developed a genetic algorithm for solutions. Caggiani et al. (2020) developed an equality-based station location model to address disparities in multi-modal mobility systems integrating bicycles and public transport. Wu et al. (2020) aimed to optimize the network of bike-sharing systems, aligning with existing bus networks.

Compared to strategic-level challenges, operational-level issues have gained increasing attention in recent years, particularly in inventory management and bike rebalancing. Inventory management focuses on distributing bike inventory across stations to enhance accessibility for users, while bike rebalancing optimizes bike transportation between stations to address imbalances between supply and demand. Regarding inventory management, Raviv and Kolka (2013) quantified service levels for users and optimized bike inventory. Datner et al. (2019) considered user utility and transit behavior across nearby stations when making inventory allocation decisions and applied a guided local search algorithm to address bike distribution planning. Schuijbroek et al. (2017) aimed to optimize station inventory levels, setting them as rebalancing targets for routing optimization, thereby integrating inventory management with bike rebalancing.

The bike rebalancing operation includes static rebalancing and dynamic rebalancing (Shui and Szeto, 2020). The static rebalancing operation focuses on the rebalancing implemented in the nighttime which the demand can be ignored. Dynamic rebalancing is typically conducted during the daytime, requiring operators to account for ongoing demand during the rebalancing process. In contrast, static bike rebalancing problems focus on optimizing rebalancing truck routes and are often formulated as vehicle routing problems (VRPs). For instance, Raviv et al. (2013) introduced two models to minimize total rebalancing costs. The cases of single and multiple rebalancing trucks were explored by Erdogan et al. (2015) and Bulhões et al. (2018), respectively. To solve VRP-based rebalancing problems, various approaches such as Benders decomposition, branch-and-cut algorithms, and heuristic methods have been applied (Erdogan et al., 2015; Zhang et al., 2017; Shui and Szeto, 2018; Bulhões et al., 2018; Lv et al., 2022). Additionally, Lv et al. (2024) incorporated both service efficiency and carbon emissions into the rebalancing process, proposing a feature correlation reinforcement clustering strategy to enhance optimization.

More related to our study, dynamic bike rebalancing is more complicated due to non-negligible demand during the rebalancing procedure. Some studies focused on optimizing the routes of rebalancing trucks during the multi-period rebalancing. For example, Chang et al. (2023) studied the rebalancing operations for free-floating systems, and they proposed a smart Predict-then-Optimize model to address the rebalancing routing problem. Moreover, Guo et al. (2024) developed a reinforcement learning approach to optimize the inventory routing problem, which considered multiple objectives like dispatch cost, waiting time, and the satisfaction of users. Besides, more research formulated the dynamic bike rebalancing problem as an inventory transshipment model, ignoring the routes of rebalancing trucks. Frade and Ribeiro (2015) incorporated two types of users, i.e., customers and subscribers in the bike-sharing systems and optimized the inventory rebalancing operation. Silva et al. (2024) presented prioritization strategies for operators to choose the stations that need to be rebalanced with priority.

2.2 Supervision in bike-sharing systems

For emerging transportation systems like ride-hailing systems and bike-sharing systems, the governments will conduct supervision strategies to improve the market environment and social welfare. For instance, Pu et al. (2020) analyzed the behaviors of the platform, passengers and drivers in online ride-hailing systems considering the supervision of the platform. Zhong et al. (2022) studied the effect of supervision strategies on ride-hailing platforms, and their results suggested that the supervision of the pricing may lower the profit of ride-hailing platforms but improve social welfare.

In bike-sharing systems, Li et al. (2021) addressed that supervision helps reduce congestion and mitigate the consumption of public resources. Wang et al. (2023) considered a supervision policy that controls the number of allocated bikes in a certain region, which can help mitigate the oversupply of bikes and address the congestion caused by crowded bikes. In addition to the supervision of oversupplying, Zhao et al. (2022) focused on the supervision of broken bikes and they employed a game framework to test various strategies of the operator and government.

2.3 Optimization approaches for uncertainty

In bike-sharing systems, operators face significant uncertainty and several methods have been developed to address it. For example, stochastic optimization, as described in Birge and Louveaux (2011), assumes that decision-makers have precise knowledge of the uncertain variables and optimizes the expected objective value based on the nominal distribution.

However, the exact probability distribution is often unknown or difficult to infer from observed data. When limited information on uncertain variables is accessible, the robust optimization (Ben-Tal and Nemirovski, 1999) and DRO (Ben-Tal et al., 2013) serve as viable alternatives. Robust optimization models uncertainty through an uncertainty set and seeks solutions that perform optimally with the worst-case scenario. Similarly, DRO constructs an ambiguity set encompassing a range of potential distributions and optimizes the expected objective function under the worst-case distribution within this set, assuming the true distribution lies within it. This method integrates empirical distributional information derived from historical data, such as moment-based details (Wiesemann et al., 2014), as well as scenario-based or probability distance-based information (Ben-Tal et al., 2013; Esfahani and Kuhn, 2018), mitigating the over-conservatism of traditional robust optimization models.

Specifically, demand uncertainty is significantly challenging for the operators, as operators must predict demand to set inventory targets and execute rebalancing operations (Yin et al., 2023). To address this issue, two-stage stochastic optimization models are commonly used, jointly optimizing bike allocation in the first stage and rebalancing operations in the second stage (Yan et al., 2017; Maggioni et al., 2019; Li et al., 2024). The robust optimization approach is also applied to mitigate demand uncertainty, providing a baseline profit for risk-averse operators. For instance, Lu (2016) formulated a robust optimization framework to address the dynamic rebalancing problem. Additionally, Chen et al. (2023) introduced a risk measure within a DRO framework to enhance the efficiency of the system. They specifically developed a lexicographic optimization model to address service equity concerns.

2.4 Research gaps

A review of existing literature reveals that, while extensive research has addressed strategic planning and operational problems, few studies have considered the influence of government supervision, which is a realistic problem for the operators. The works of Zhao et al. (2022) and Wang et al. (2023) explored the effects of government regulation on bike-sharing systems. However, these studies only focused on analyzing the impact of supervision rather than the operation strategy of the bike-sharing systems when facing supervision. Furthermore, the uncertainty of government supervision has not been addressed in current research.

In this study, we jointly optimize the strategic and operational decisions in bike-sharing systems, incorporating the uncertain government supervision strategy. Specifically, the supervision for oversupplying of bikes and malfunctioning bikes is included. To address the demand and the government supervision strategy uncertainty, we propose a DRO model that leverages both moment and sample information. To the best of our knowledge, this is the first study to optimize the strategic and operational decisions of bike-sharing systems under government supervision. Moreover, we derive managerial insights for both bike-sharing operators and the government supervisors based on numerical results.

3 Problem description and model formulation

In this paper, we jointly optimize the depot location and the daytime operation in free-floating bike-sharing systems under the supervision of the local government. We assume that the whole operational area of the free-floating bike-sharing system has been divided into smaller operational regions, and the set of operational regions is defined as \mathcal{R}=\{1,2,\dots ,R\}. We denote \mathcal{T} as the set of daily operational periods, and we have \mathcal{T}=\{1,2,\dots ,T\}. In each period t, travel demand may occur between any origin-destination (O-D) pair \left(i,j\right) for any i,j\in \mathcal{R}, and we denote travel demand from region i to region j during period t as D_{i,j,t}.

The operator should allocate bikes among all operational regions to satisfy the demand. We denote W as the number of bikes allocated across the whole system at the beginning of one day, and use n_{i,t}^{\mathrm{a}} to represent available bikes in region i at the beginning of period t. Specifically, when t=1, n_{i,1}^{\mathrm{a}} represents initial bike allocations in region i at the beginning of the day. The bikes allocated among these operational regions can be used by users to satisfy their travel demand. We use d_{i,j,t} to represent satisfied demand from region i to region j during period t by users, and the unit revenue of satisfying a travel demand is r. Considering the limited supply of bikes in the systems, not all travel demand may be satisfied. We define u_{i,j,t} as unmet demand from region i to region j in period t. The unmet demand can cause dissatisfaction among users and even damage the reputation of the operator, affecting the long-term profit of the system. Thus, each unsatisfied demand of the user will incur a penalty cost of c_1.

To avoid improper operations, the government may implement supervision of the bike-sharing systems. We consider two types of supervision in this paper, i.e., supervision of malfunctioning bikes and oversupplying bikes, which have been employed in many cities with bike-sharing systems (Shenzhen Municipal Transportation Bureau, 2017). The two supervision strategies will be introduced in detail as follows:

Supervision of malfunctioning bikes. To keep the bike-sharing systems maintaining operational efficiency, the government may conduct supervision of bike-sharing systems and the supervision of malfunctioning bikes is one of the key measures because the malfunction of bikes is a significant operational challenge. Bike malfunctions typically occur after usage, and we assume that the malfunction probability of available bikes after usage is \alpha (Jin et al., 2022). We use n_{i,t}^{\mathrm{b}} to denote malfunctioning bikes in region i at the beginning of period t. Malfunctioning bikes cannot satisfy demand and may lead to user dissatisfaction, as riders might walk a long distance only to find an unusable bike. Thus, users can report the malfunctioning bikes to the government supervisor. We assume that in each period t, each malfunctioning bike will be punished with a penalty cost c_3 by the supervisor.

Supervision of oversupplying bikes. In addition to the supervision of malfunctioning bikes, the supervision of oversupplying bikes is also frequently implemented to restrict the supply of bikes to a reasonable level to avoid the congestion caused by oversupplying bikes. We denote h_{i,t} as the supervision strategy of the supervisor, i.e., h_{i,t}=1 means that the supervisor chooses to supervise region i in period t and otherwise h_{i,t}=0.

To avoid congestion caused by the oversupply of bikes, the supervisor usually stipulates a maximum bike inventory level U_{i,t} for each region i\in \mathcal{R} and period t\in \mathcal{T}. Note that the unavailable malfunctioning bikes will also account for the inventory of bikes in each region, thus the maximum inventory level U_{i,t} includes both available bikes and malfunctioning bikes. We assume that the supervisor will arrive at region i and monitor the number of bikes in this region at the beginning of period t if h_{i,t}=1. If the actual number of bikes in region i at the beginning of period t is larger than the maximum inventory level U_{i,t}, i.e., n_{i,t}^{\mathrm{a}}+n_{i,t}^{\mathrm{b}}\ge U_{i,t}, the operator of the bike-sharing system will be punished with a penalty cost for each bike exceeding the maximum inventory level. We define the number of bikes exceeding the maximum inventory level in region i at the beginning of period t as e_{i,t}, i.e., e_{i,t}=\mathrm{m}\mathrm{a}\mathrm{x}{\{n}_{i,t}^{\mathrm{a}}+n_{i,t}^{\mathrm{b}}-U_{i,t},0\}, and the penalty cost for each excess allocated bike is c_2. Note that for any region i\in \mathcal{R} and period t\in \mathcal{T}, the penalty cost for oversupplying supervision in this region and period only occurs if the supervisor chooses to supervise region i at the beginning of period t, i.e., h_{i,t}=1. And if h_{i,t}=0, the existence of the oversupplying bikes in region i and period t will not be punished because the supervisor cannot discover these improper operations. We define the number of bikes exceeding the maximum inventory level that is punished by the supervisor in region i at the beginning of period t as v_{i,t}.

Depot location and bike replenishment. To avoid the punishment of supervised malfunctioning and oversupplied bikes while enhancing service quality, bike-sharing operators typically establish multiple depots within operational regions to store additional available bikes. During the daytime operation, the operator can move malfunctioning bikes from operational regions to the depots, and replenish available bikes from the depots to operational regions to maintain a high proportion of available bikes and avoid the punishment of malfunctioning bikes. At night, the malfunctioning bikes that have been collected in the depots during the daytime can be repaired or replaced with available bikes by the operator, so that these available bikes in the depots can be used to replenish malfunctioning bikes in the next day. Suppose that the operator can construct P homogeneous depots, and the capacity of available and malfunctioning bikes in each depot are Q^{\mathrm{a}} and Q^{\mathrm{b}}. Thus, the bikes can be allocated among operational regions and depots. The operator should determine the location of these depots to reduce the cost of transporting bikes. We use binary variable x_i to represent the depot location decision. If a depot is located in region i\in \mathcal{R} then we have x_i=1, and otherwise x_i=0. For brevity, we use depot i to represent the depot located in region i. In addition, we define y_{i,j,t} as available bikes moved from depot i to region j during period t. Similarly, we define {\mathit{\text{z}}}_{i,j,t} as malfunctioning bikes moved from region i to depot j. Bike transportation between regions and depots is implemented by the truck fleet of the operator with a capacity of C. The cost of transferring a bike from a region (depot) i to depot (region) j is c_{4,i,j}, which corresponds to the distance between region i and j. The inventory of available bikes and malfunctioning bikes stored in depot i at period t is represented as q_{i,t}^{\mathrm{a}} and q_{i,t}^{\mathrm{b}} respectively.

Bike rebalancing. In addition to bike replenishment from depots, the operator usually implements bike rebalancing operations to address the imbalance of bike supply and avoid the punishment of oversupplied bikes in some regions. We suppose that the operator conducts rebalancing operations at the beginning of each period t, and we denote k_{i,j,t} as available bikes rebalanced from region i to region j. The bike rebalancing operation is also implemented by the truck fleet, thus the rebalancing cost of transferring one bike from region i to region j is also c_{4,i,j}.

The detailed description of the notations used in this paper is presented in Table 1.

In the objective Eq. (1a), we optimize the profits in system, which is reflected by the revenue minus total cost. The revenue is derived from satisfying the demand of users, and the cost includes the penalty cost of unsatisfied demand, oversupplying and malfunctioning bikes. Constraints (1b)–(1e) are the depot location and initial bike allocation constraints. Constraint (1b) means that there are P depots to be located. Constraints (1c)–(1e) mean that there are a total number of W available bikes allocated among all regions and depots at the beginning of one day, i,e., t=1, and no malfunctioning bikes exist in the system at the beginning of one day. Constraints (1f) and (1g) mean that there is no inventory in depot i if no depot is located in region i, and the inventory of both available and malfunctioning bikes in depot i cannot exceed the corresponding capacity Q^{\mathrm{a}} and Q^{\mathrm{b}}.

(1a)\begin{array}{c}\begin{array}{c}\mathrm{m}\mathrm{a}\mathrm{x}{\displaystyle \sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}r{d}_{i,j,t}-\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{c}_1{u}_{i,j,t}-\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{c}_2{v}_{i,t}}\\ -{\displaystyle \sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{c}_3{n}_{i,t}^{\mathrm{b}}-\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{c}_{4,i,j}\left(k_{i,j,t}+y_{i,j,t}+{\mathit{\text{z}}}_{i,j,t}\right),}\end{array}\end{array}

(1b)\begin{array}{c}\text{s.t.} \sum _{i\in \mathcal{R}}{x}_i=P,\end{array}

(1c)\begin{array}{c}\sum _{i\in \mathcal{R}}{n}_{i,1}^{\mathrm{a}}+\sum _{i\in \mathcal{R}}{q}_{i,1}^{\mathrm{a}}=W,\end{array}

(1d)\begin{array}{c}{n}_{i,1}^{\mathrm{b}}=0,\mathrm{\forall }i\in \mathcal{R},\end{array}

(1e)\begin{array}{c}{q}_{i,1}^{\mathrm{b}}=0,\mathrm{\forall }i\in \mathcal{R},\end{array}

(1f)\begin{array}{c}{q}_{i,t}^{\mathrm{a}}\le Q^{\mathrm{a}}{x}_i,\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},\end{array}

(1g)\begin{array}{c}{q}_{i,t}^{\mathrm{b}}\le Q^{\mathrm{b}}{x}_i,\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T}.\end{array}

Constraints (2a)–(2h) are the operational constraints in the system, including demand satisfaction, bike rebalancing, malfunctioning bike transportation and available bike replenishment between depots and regions. Specifically, constraint (2a) ensures that in each period, the bikes taken by users and for rebalancing from region i should be less than the available bike inventory at the beginning of that period. Constraint (2b) limits the malfunctioning bikes collected from region i to depots at period t to malfunctioning bikes in that region at the same time. Constraint (2c) ensures that the available bikes replenished from depot i do not exceed the depot’s available bike inventory. Constraints (2d) and (2e) mean that malfunctioning bikes can be collected back to depot i and available bikes can be replenished from depot i only if a depot is located in region i. Constraint (2f) restricts the number of satisfied demand should not exceed the total demand for any region i,j and period t, and Constraint (2g) is the restriction of unsatisfied demand. In addition, Constraint (2h) is the capacity constraint of the rebalancing truck fleet. Note that the malfunctioning bike collection, available bike replenishment, and bike rebalancing are all implemented by the truck fleet. Accordingly, these operations in each period are restricted by the capacity C of the truck fleet.

(2a)\begin{array}{c}\sum _{j\in \mathcal{R}}\left(d_{i,j,t}+k_{i,j,t}\right)\le n_{i,t}^{\mathrm{a}},\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},\end{array}

(2b)\begin{array}{c}\sum _{j\in \mathcal{R}}{\mathit{\text{z}}}_{i,j,t}\le n_{i,t}^{\mathrm{b}},\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},\end{array}

(2c)\begin{array}{c}\sum _{j\in \mathcal{R}}{y}_{i,j,t}\le q_{i,t}^{\mathrm{a}},\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},\end{array}

(2d)\begin{array}{c}{y}_{i,j,t}\le M{x}_i,\mathrm{\forall }i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},\end{array}

(2e)\begin{array}{c}{\mathit{\text{z}}}_{j,i,t}\le M{x}_i,\mathrm{\forall }i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},\end{array}

(2f)\begin{array}{c}{d}_{i,j,t}\le D_{i,j,t},\mathrm{\forall }i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},\end{array}

(2g)\begin{array}{c}{u}_{i,j,t}\ge D_{i,j,t}-d_{i,j,t},\mathrm{\forall }i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},\end{array}

(2h)\begin{array}{c}\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\left(k_{i,j,t}+y_{i,j,t}+{\mathit{\text{z}}}_{i,j,t}\right)\le C,\mathrm{\forall }t\in \mathcal{T}.\end{array}

Constraints (3a) and (3b) are related to the supervision of oversupplied bikes. Constraint (3a) gives oversupplied bikes in region i at period t. Constraint (3b) gives oversupplied bikes that are punished by the supervisor in region i at period t.

(3a)\begin{array}{c}\begin{array}{c}\end{array}{e}_{i,t}\ge n_{i,t}^{\mathrm{a}}+n_{i,t}^{\mathrm{b}}-U_{i,t},\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},\end{array}

(3b)\begin{array}{c}{v}_{i,t}\ge h_{i,t}{e}_{i,t},\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T}.\end{array}

Constraints (4a)–(4d) illustrate the inventory variation constraints between two consecutive periods for any operational regions and depots. Specifically, Constraints (4a) and (4b) illustrate the inventory variation of depots. Constraint (4a) means that the available bikes in depot i at period t equals the available bikes in depot i at period t-1 minus bikes replenished to all regions from depot i in period t-1. And constraint (4b) stipulates the inventory variation of malfunctioning bikes between period t-1 and t, i.e., for any depot i, the inventory of malfunctioning bikes at period t equals the inventory at period t-1 plus the new collected malfunctioning bikes from all regions in period t-1. Furthermore, Constraints (4c) and (4d) illustrate inventory variation of operational regions. Constraint (4c) gives the inventory of available bikes in region i at period t. Compared with n_{i,t}^{\mathrm{a}}, the outflow of available bikes to region i includes the bikes used to satisfy demand and rebalanced to other regions from region i. The inflow includes the bikes rebalanced from other regions, replenished from depots and the non-malfunctioning proportion (1-\alpha) of bikes after usage ridden by users from other regions to region i. Similarly, Constraint (4d) shows the inventory variation of malfunctioning bikes in region i between period t-1 and t. The outflow of malfunctioning bikes is bikes that are collected to depots and the inflow is the number of malfunctioning bikes after the cycling of users from other zones. Constraints (4e) and (4f) are the domains of all variables in the model.

(4a)\begin{array}{c}{q}_{i,t}^{\mathrm{a}}=q_{i,t-1}^{\mathrm{a}}-\sum _{j\in \mathcal{R}}{y}_{i,j,t-1},\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},t\ge 2,\end{array}

(4b)\begin{array}{c}{q}_{i,t}^{\mathrm{b}}=q_{i,t-1}^{\mathrm{b}}+\sum _{j\in \mathcal{R}}{\mathit{\text{z}}}_{j,i,t-1,}\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},t\ge 2,\end{array}

(4c)\begin{array}{c}\begin{array}{rl}{n}_{i,t}^{\mathrm{a}}=& n_{i,t-1}^{\mathrm{a}}-\sum _{j\in \mathcal{R}}\left(d_{i,j,t-1}+k_{i,j,t-1}\right)\\ & +\sum _{j\in \mathcal{R}}\left(\left(1-\alpha \right)d_{j,i,t-1}+k_{j,i,t-1}+y_{j,i,t-1}\right),\\ & \mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},t\ge 2,\end{array}\end{array}

(4d)\begin{array}{c}{n}_{i,t}^{\mathrm{b}}=n_{i,t-1}^{\mathrm{b}}-\sum _{j\in \mathcal{R}}{\mathit{\text{z}}}_{i,j,t-1}+\alpha \sum _{j\in \mathcal{R}}{d}_{j,j,t-1},\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},t\ge 2,\end{array}

(4e)\begin{array}{c}{x}_i\in \left\{0,1\right\},\mathrm{\forall }i\in \mathcal{R},\end{array}

(4f)\begin{array}{c}\begin{array}{c}{q}_{i,t}^{\mathrm{a}},q_{i,t}^{\mathrm{b}},n_{i,t}^{\mathrm{a}},n_{i,t}^{\mathrm{b}},e_{i,t},v_{i,t},k_{i,j,t},d_{i,j,t},{\mathit{\text{z}}}_{i,j,t},y_{i,j,t},u_{i,j,t}\in {\mathbb{Z}}^{+},\\ \mathrm{\forall }i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T}.\end{array}\end{array}

4 Optimization models considering uncertainty of supervision strategy and demand

Though the above model provides a tool for the operators in bike-sharing systems under supervision, it relies on the accurate estimation of both the travel demand of users and the supervision strategy of the supervisor. However, in practice, the travel demand may never be estimated accurately, and the operators also have no way to obtain the supervision strategy of the supervisor. Using inaccurately estimated parameters as input to the deterministic model will lead to poor out-of-sample performance. Thus, the operators should also consider the uncertainty of both travel demand and supervision strategy in bike-sharing systems when making operational decisions.

To address the concern of uncertainty, in Subsection 4.1, we introduce a stochastic optimization model that can benefit from historical samples. Furthermore, we propose an ambiguity set with both moment-based and sample-based information and introduce a two-stage DRO model in Subsection 4.2.

4.1 Stochastic optimization model

To tackle the uncertainty of supervision strategy and demand uncertainty, the stochastic optimization approach is commonly used to optimize the expected objective based on a given probability distribution of the uncertain variable. We denote \mathbb{P} as the joint probability distribution of supervision strategy by the supervisor (\tilde{\mathbf{h}}) and the travel demand from users (\tilde{\mathbf{D}}). Note that the true distribution \mathbb{P} is not easily obtained by the operator in practice. As an alternative, the true distribution \mathbb{P} can be approximated by the empirical distribution \hat{\mathbb{P}}, which can be estimated according to the historical operational data in the bike-sharing system. Suppose that the operator of the bike-sharing system is accessible to a series of scenarios corresponding to the observation of historical data, and the set of scenarios is denoted as \mathcal{S}=\{1,2,\dots ,S\}. The empirical distribution \hat{\mathbb{P}} can be represented as \hat{\mathbb{P}}\left[\left(\tilde{\mathbf{h}},\tilde{\mathbf{D}}\right)=\left({\tilde{\mathbf{h}}}_s,{\tilde{\mathbf{D}}}_s\right)\right]=\frac{1}{S}, in which \left({\tilde{\mathbf{h}}}_s,{\tilde{\mathbf{D}}}_s\right) is a realized historical data under scenario s\in \mathcal{S}.

In the two-stage stochastic optimization model, the depot location and initial bike allocation decisions are optimized in the first stage. The second-stage decisions are operational-level decisions, including other operational decisions in the deterministic model. For the purpose of modeling, we additionally define decision variables N{n}_i^{\mathrm{a}} and N{q}_i^{\mathrm{a}} as the initial allocation of available bikes at the beginning of each day in region i and depot i respectively, i.e., we have N{n}_i^{\mathrm{a}}=n_{i,1}^{\mathrm{a}} and N{q}_i^{\mathrm{a}}=q_{i,1}^{\mathrm{a}} for any i\in \mathcal{R}. The two-stage stochastic optimization Model (5) is formulated as follows:

(5a)\begin{array}{c}\mathrm{m}\mathrm{a}\mathrm{x}{\mathbb{E}}_{\hat{\mathbb{P}}}\left[H\left(\mathbf{x},{\mathbf{N}\mathbf{n}}^{\mathrm{a}},{\mathbf{N}\mathbf{q}}^{\mathrm{a}},\tilde{\mathbf{h}},\tilde{\mathbf{D}}\right)\right],\end{array}

(5b)\begin{array}{c}\text{s.t.} \sum _{i\in \mathcal{R}}{x}_i=P,\end{array}

(5c)\begin{array}{c}\sum _{i\in \mathcal{R}}N{n}_i^{\mathrm{a}}+\sum _{i\in \mathcal{R}}N{q}_i^{\mathrm{a}}=W,\end{array}

(5d)\begin{array}{c}N{q}_i^{\mathrm{a}}\le M{x}_i,\mathrm{\forall }i\in \mathcal{R},\end{array}

(5e)\begin{array}{c}{x}_i\in \left\{0,1\right\},\mathrm{\forall }i\in \mathcal{R},\end{array}

(5f)\begin{array}{c}N{n}_i^{\mathrm{a}},N{q}_i^{\mathrm{a}}\in {\mathbb{Z}}^{+},\mathrm{\forall }i\in \mathcal{R}.\end{array}

For simplicity, we define vector \nu to represent all first-stage variables including the depot location and bike allocation decisions, i.e., \mathit{\nu }=(x_i,N{n}_i^{\mathrm{a}},N{q}_i^{\mathrm{a}}{)}_{i\in \mathcal{R}}, and replace H\left(\mathbf{x},{\mathbf{N}\mathbf{n}}^{\mathrm{a}},{\mathbf{N}\mathbf{q}}^{\mathrm{a}},\tilde{\mathbf{h}},\tilde{\mathbf{D}}\right) with H\left(\mathit{\nu },\tilde{\mathbf{h}},\tilde{\mathbf{D}}\right). The feasible set of \nu is defined as \mathcal{X} for convenience, which includes Constraints (5b)–(5f). We present the second-stage problem H\left(\mathit{\nu },\tilde{\mathbf{h}},\tilde{\mathbf{D}}\right) as follows:

(6a)\begin{array}{c}\begin{array}{rl}H\left(\nu ,\tilde{\mathit{h}},\tilde{\mathit{D}}\right)=& \mathrm{m}\mathrm{a}\mathrm{x}\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}r{d}_{i,j,t}-\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{c}_1{u}_{i,j,t}\\ & -\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{c}_2{v}_{i,t}-\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{c}_3{n}_{i,t}^{\mathrm{b}}-\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}\\ & c_{4,i,j}(k_{i,j,t}+y_{i,j,t}+z_{i,j,t}),\end{array}\end{array}

\text{s.t.}\text{Constraints (1}\text{d}\text{) - (1}\text{g}\text{), (2}\text{a}\text{) - (2}\text{e}\text{), (2}\text{h}\text{), (3}\text{a}\text{), (4}\text{a}\text{) - (4}\text{d}\text{)}\text{,}

(6b)\begin{array}{c}{d}_{i,j,t}\le {\tilde{D}}_{i,j,t},\mathrm{\forall }i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},\end{array}

(6c)\begin{array}{c}{u}_{i,j,t}\ge {\tilde{D}}_{i,j,t}-d_{i,j,t},\mathrm{\forall }i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},\end{array}

(6d)\begin{array}{c}{v}_{i,t}\ge {\tilde{h}}_{i,t}{e}_{i,t,}\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},\end{array}

(6e)\begin{array}{c}N{n}_i^{\mathrm{a}}=n_{i,1}^{\mathrm{a}},\mathrm{\forall }i\in \mathcal{R},\end{array}

(6f)\begin{array}{c}N{q}_i^{\mathrm{a}}=q_{i,1}^{\mathrm{a}},\mathrm{\forall }i\in \mathcal{R}.\end{array}

Model (5) optimizes the expected profits given the empirical distribution. Suppose that each scenario s\in \mathcal{S} has the same probability, i.e., p_s={\displaystyle \frac{1}{S}} for any s\in \mathcal{S}, the stochastic optimization model is reformulated into Model (7), where H^s\left(\nu ,{\mathbf{h}}^{\mathbf{s}},{\mathbf{D}}^{\mathbf{s}}\right) represents the second-stage problem for each scenario s\in \mathcal{S}.

(7a)\begin{array}{c}\underset{\nu }{\mathrm{m}\mathrm{a}\mathrm{x}}\sum _{s\in \mathcal{S}}{p}_s\left[H^s\left(\mathit{\nu },{\mathbf{h}}^{\mathbf{s}},{\mathbf{D}}^{\mathbf{s}}\right)\right],\end{array}

(7b)\begin{array}{c}\text{s.t.}\mathit{\nu }\in \mathcal{X}.\end{array}

Specifically, the constraints that are related to the uncertain variables in problem H^s\left(\nu ,{\mathbf{h}}^{\mathbf{s}},{\mathbf{D}}^{\mathbf{s}}\right) are represented as the scenario-wise constraints set (8), where D_{i,j,t}^s and h_{i,t}^s are denoted as the related demand and supervision strategy in historical sample s\in \mathcal{S}.

(8a)\begin{array}{c}{d}_{i,j,t}^s\le D_{i,j,t}^s,\mathrm{\forall }i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},s\in \mathcal{S},\end{array}

(8b)\begin{array}{c}{u}_{i,j,t}^s\ge D_{i,j,t}^s-d_{i,j,t}^s,\mathrm{\forall }i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},s\in \mathcal{S},\end{array}

(8c)\begin{array}{c}{v}_{i,t}^s\ge h_{i,t}^s{e}_{i,t,}^s\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},s\in \mathcal{S}.\end{array}

Though employing the historical samples, the stochastic optimization approach also requires sufficient data of the uncertain demand to estimate the true distribution (He et al., 2020). As a result, the stochastic optimization model using empirical distribution may still yield unreliable solutions with poor performance when the empirical distribution deviates from the true distribution (Long et al., 2023). In practice, the operator cannot obtain the full information about the supervision strategy and travel demand. In most cases, only finite samples are observed, and only limited information about the true distribution is accessible, which is not sufficient to estimate the true distribution. Smith and Winkler (2006) discovered that the optimal objectives and solutions derived by the empirical optimization are unreliable, referring to the “optimizer’s curse”. To tackle this issue, the operator may seek an approach that can benefit from partial historical information and derive reliable solutions against any possible distributions. Thus, we propose a DRO framework that can benefit from finite historical samples and partial distributional information to address these concerns in the following sections.

4.2 DRO model

Different from the stochastic optimization model that requires perfect knowledge of the uncertain variables, the DRO model assumes that the distribution of uncertain variables is unknown to the decision-maker and contained in an ambiguity set. The DRO model aims to optimize the expected objective with the worst-case distribution in the ambiguity set to avoid risk.

Before formulating the DRO model, we first propose the ambiguity set. We first define the joint distribution of {\tilde{D}}_{i,j,t} and {\tilde{h}}_{i,t} as \mathbb{P}, and it is assumed to be contained within an ambiguity set \mathbb{F}. The operators of bike-sharing systems may have the daily operational demand data and the observations of the supervision strategy. To better benefit from the historical data and avoid operational risk, we formulate an ambiguity set that includes both moment-based and sample-based information as follows:

\mathbb{F}=\left\{\mathbb{P}\in {\mathcal{P}}_0|\begin{array}{lll}& {\displaystyle \sum _{s\in \mathcal{S}}{\tilde{p}}_s=1,}& \\ & \left|{\mathbb{E}}_{\mathbb{P}}\left({\tilde{h}}_{i,t}\right)-{\mu }_{i,t}\right|\le {ϵ}_{i,t},& \mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},\\ & \left|{\mathbb{E}}_{\mathbb{P}}\left({\tilde{D}}_{i,j,t}\right)-{\theta }_{i,j,t}\right|\le {\sigma }_{i,j,t},& \mathrm{\forall }i,j\in \mathcal{R},t\in \mathcal{T},\\ & {\tilde{p}}_s\in \left[{\underset{\_}{p}}_s,{\overline{p}}_s\right],& \mathrm{\forall }s\in \mathcal{S}.\end{array}\right\}

In ambiguity set \mathbb{F}, {\tilde{p}}_s represents the probability of scenario s under probability distribution \mathbb{P}. In addition, {\mu }_{i,t} and {ϵ}_{i,t} are the mean and mean absolute deviation of the supervision frequency in region i at period t under the empirical distribution, respectively. Similarly, {\theta }_{i,j,t} and {\sigma }_{i,j,t} are the mean and mean absolute deviation of uncertain travel demand of users from region i to region j during period t under the empirical distribution, respectively. Thus, we have {\mathbb{E}}_{\mathbb{P}}\left({\tilde{h}}_{i,t}\right)=\sum _{s\in \mathcal{S}}{\tilde{p}}_s{h}_{i,t}^s and {\mathbb{E}}_{\mathbb{P}}\left({\tilde{D}}_{i,j,t}\right)= \sum _{s\in \mathcal{S}}{\tilde{p}}_s{D}_{i,j,t}^s. The second and third constraints represent the moment-based information of any distribution \mathbb{P}. These constraints restrict the feasible region of the probability distributions \mathbb{P} in the ambiguity set \mathbb{F}.

Notice that the supervision strategy h_{i,t} remains as binary uncertain variable in the ambiguity set. Because the mathematical expectation {\mathbb{E}}_{\mathbb{P}}\left({\tilde{h}}_{i,t}\right) is represented as \sum _{s\in \mathcal{S}}{\tilde{p}}_s{h}_{i,t}^s, which is related to the historical observations of the binary supervision strategy h_{i,t}^s for each s\in \mathcal{S}. Moreover, {\tilde{p}}_s is also restricted within an interval \left[{\underset{\_}{p}}_s,{\overline{p}}_s\right], which represent the sample-based information of \mathbb{P}. Specifically, the parameters related to the empirical distribution \hat{\mathbb{P}}, i.e., the mean and mean absolute deviation of the supervision frequency and travel demand, are obtained from historical samples, which can be calculated as follows:

(9a)\begin{array}{c}{\mu }_{i,t}=\frac{1}{S}\sum _{s\in \mathcal{S}}{h}_{i,t}^s,\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},\end{array}

(9b)\begin{array}{c}{ϵ}_{i,t}=\frac{1}{S}\sum _{s\in \mathcal{S}}|h_{i,t}^s-{\mu }_{i,t}|,\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},\end{array}

(9c)\begin{array}{c}{\theta }_{i,j,t}=\frac{1}{S}\sum _{s\in \mathcal{S}}{D}_{i,j,t}^s,\mathrm{\forall }i,j\in \mathcal{R},t\in \mathcal{T},\end{array}

(9d)\begin{array}{c}{\sigma }_{i,j,t}=\frac{1}{S}\sum _{s\in \mathcal{S}}|D_{i,j,t}^s-{\theta }_{i,j,t}|,\mathrm{\forall }i,j\in \mathcal{R},t\in \mathcal{T}.\end{array}

In Eq. (9), h_{i,t}^s and D_{i,j,t}^s represent the related supervision strategy and travel demand in historical sample s\in \mathcal{S}. In addition, the decision-maker can calibrate the ranging parameters {\underset{\_}{p}}_s and {\overline{p}}_s through such pre-experiments to choose proper parameters. And we test the effect of ranging parameters on the model performance in Subsection 5.2. Compared with the stochastic optimization Model (5) that only relies on the empirical distribution, the DRO model additionally incorporates the moment distributional information and allows the probability of each scenario to vary in a certain range. Incorporating the additional moment-based information, the DRO model can become more risk-averse compared with the stochastic optimization model.

The bike-sharing operation model considering supervision and demand uncertainty is presented as follows, given the ambiguity set \mathbb{F}:

(10a)\begin{array}{c}\underset{\nu }{\mathrm{m}\mathrm{a}\mathrm{x}}\underset{\mathbb{P}\in \mathbb{F}}{\mathrm{i}\mathrm{n}\mathrm{f}}{\mathbb{E}}_{\mathbb{P}}\left[H\left(\mathit{\nu },\tilde{\mathbf{h}},\tilde{\mathbf{D}}\right)\right],\end{array}

(10b)\begin{array}{c}\text{s.t.}\mathit{\nu }\in \mathcal{X}.\end{array}

In Model (10), we optimize the profit of the bike-sharing system over the worst-case distribution \mathbb{P} in ambiguity set. However, the two-stage DRO Model (10) is non-convex due to ambiguity set \mathbb{F} and the inner minimization problem H\left(\mathit{\nu },\tilde{\mathbf{h}},\tilde{\mathbf{D}}\right). In Proposition 1, we show the reformulation of Model (10).

Proposition 1. The DRO Model (10) can be reformulated into the following formulation according to the ambiguity set \mathbb{F}:

(11a)\begin{array}{c}\begin{array}{rl}& \mathrm{m}\mathrm{a}\mathrm{x}\lambda -\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}\left({\mu }_{i,t}+{ϵ}_{i,t}\right){\eta }_{i,t}^1+\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}\left({\mu }_{i,t}-{ϵ}_{i,t}\right){\eta }_{i,t}^2\\ & -\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}\left({\theta }_{i,j,t}+{\sigma }_{i,j,t}\right){\varphi }_{i,j,t}^1\\ & +\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}\left({\theta }_{i,j,t}-{\sigma }_{i,j,t}\right){\varphi }_{i,j,t}^2+\sum _{s\in \mathcal{S}}\left({\underset{\_}{p}}_s{\psi }_s^1-{\overline{p}}_s{\psi }_s^2\right),\end{array}\end{array}

(11b)\begin{array}{c}\begin{array}{rl}& \text{s.t.}\lambda -\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{h}_{i,t}^s{\eta }_{i,t}^1+\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{h}_{i,t}^s{\eta }_{i,t}^2-\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{D}_{i,j,t}^s{\varphi }_{i,j,t}^1\\ & +\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{D}_{i,j,t}^s{\varphi }_{i,j,t}^2+{\psi }_s^1-{\psi }_s^2\le H\left(\mathit{\nu },{\mathbf{h}}^{\mathbf{s}},{\mathbf{D}}^{\mathbf{s}}\right),\\ & \mathrm{\forall }s\in \mathcal{S},\end{array}\end{array}

(11c)\begin{array}{c}{\eta }_{i,t}^1,{\eta }_{i,t}^2,{\varphi }_{i,j,t}^1,{\varphi }_{i,j,t}^2,{\psi }_s^1,{\psi }_s^2\ge 0,\mathrm{\forall }s\in \mathcal{S},i,j\in \mathcal{R},t\in \mathcal{T},\end{array}

(11d)\begin{array}{c}\mathit{\nu }\in \mathcal{X}.\end{array}

The proof of all propositions is shown in Appendixes. We note that the proposed two-stage DRO model can be reformulated into an equivalent MIP model through Proposition 1 and solved by commercial solvers directly. In addition, we will illustrate that the stochastic optimization model proposed in Section 4.1 can be regarded as a special case of the DRO model in Proposition 2.

Proposition 2. Suppose that {\overline{p}}_s={\underset{\_}{p}}_s={\displaystyle \frac{1}{S}} for any s\in \mathcal{S} and {ϵ}_{i,t}={\sigma }_{i,j,t}=0 for any i,j\in \mathcal{R} and t\in \mathcal{T} in the ambiguity set \mathbb{F}, then \mathbb{F} will shrink to a singleton that only includes the empirical distribution \hat{\mathbb{P}} and the DRO model will reduce to the stochastic optimization model proposed in Section 4.1.

Thus, our proposed DRO model can be regarded as a risk-averse extension of the stochastic optimization model by incorporating additional moment-based information and allowing the probability of each scenario s\in \mathcal{S} to vary in a certain range. A wider range of \left[{\underset{\_}{p}}_s,{\overline{p}}_s\right] indicates that the ambiguity set can include more possible distributions and the DRO model will become more conservative. The detailed formulation of the reformulated equivalent deterministic model of the DRO model is shown in the Appendix.

5 Numerical experiment

In this section, we conduct numerical experiments to test the performance of the proposed model and aim to derive some managerial insights for the operators. Specifically, in Section 5.1, we describe the data used in the numerical experiments and the parameter setting. In Section 5.2, we assess the out-of-sample performance of the proposed DRO model. In Section 5.3, we conduct the effect analysis of the model and provide managerial insights for the operators in bike-sharing systems with supervision. The numerical experiments are implemented in Python 3.7 using Gurobi 9.0.3 with a Windows 10 PC equipped with an Intel Core i7-10750H CPU (2.60 GHz) and 16 GB of RAM.

5.1 Data description and parameter setting

Data description: In this study, we utilize travel data from Mobike, one of the largest free-floating bike-sharing systems in China, collected in Beijing, and a high-demand urban area was selected for analysis. The selected area is partitioned into 36 regions, i.e., R=36, and each region is about 1 km². The travel records are mapped to these regions to construct the O-D demand network. Specifically, trip data from May 10 to May 16, 2017, is used to generate scenario-based demand in our model, yielding S=7. The data reveals substantial demand variability, with daily trip volumes ranging from a minimum of 5,808 to a maximum of 8,534. Fig. 2 illustrates the hourly distribution of average daily demand, showing that over 85% of trips occur between 6:00 AM and 8:00 PM. Based on this observation, we define a 14-hour operational window for the numerical study. Furthermore, the data indicates that most trips are completed within two hours. Therefore, we divide the operational period into 2-hour intervals, resulting in T=7 time periods.

For the government supervision strategy h_{i,t}, we assume that the government supervisor will randomly select S_t regions to supervise in each period t\in \mathcal{T}. Specifically, a time-varying adaptive supervision strategy is employed in the numerical experiment to generate the supervision strategy h_{i,t}, i.e., the supervision intensity in each period is related to the demand in this period, and we set S_t\in \{3,5,4,5,4,5,4\} for each period t\in \mathcal{T}. A higher historical travel demand during period t results in higher supervision intensity S_t in this period. Moreover, regions with higher travel demand are more likely to be supervised, which reflects the spatial and temporal adaptation of the supervision strategy. In addition, the maximum inventory level in each region stipulated by the supervisor U_{i,t} is set according to the travel demand density in this region. We first set the basic setting U_{i,t}=U_{i,t}^{\mathrm{b}\mathrm{a}\mathrm{s}\mathrm{i}\mathrm{c}}, in which U_{i,t}^{\mathrm{b}\mathrm{a}\mathrm{s}\mathrm{i}\mathrm{c}} varies from 50 to 70, related to the demand density of each region and period. In Section 5.3, we also test the effect of the maximum inventory level U_{i,t} on the model performance by setting U_{i,t}=U_{\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{e}\mathrm{n}\mathrm{t}}\times U_{i,t}^{\mathrm{b}\mathrm{a}\mathrm{s}\mathrm{i}\mathrm{c}}, where U_{\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{e}\mathrm{n}\mathrm{t}} represents scale parameter.

Parameter settings: The budget of depots P is 3, and for each depot, the capacity for bikes Q^{\mathrm{a}} and Q^{\mathrm{b}} is set to 100. The capacity of the rebalancing fleet C is 50. In addition, the unit revenue of satisfying the demand of a user r is 2. The penalty cost of unmet demand c_1, the penalty cost of oversupplied bikes c_2, and the penalty cost of malfunctioning bikes c_3 are set as 5, 10, and 5, respectively. The rebalancing cost c_{4,i,j} corresponds to the travel distance between regions i and j, equaling 0.5 times the Manhattan distance between the two regions. The available bike’s malfunction probability \alpha after it finishes a trip is 0.01. We focus on the performance of the proposed models under various budgets of available bikes W. Thus, we set W\in \{1300,1400,1500, 1600,1700\}.

Simulation settings: To evaluate the performance of the proposed DRO model, we use the deterministic model (DM) as a benchmark. Specifically, the demand and supervision strategy parameters in the DM are set as the average value of the historical data. The detailed formulations of the DM and DRO models are shown in Section 3 and Appendix. In addition, the sample average approximation (SAA) model and robust optimization (RO) models with different uncertain sets are selected as testing models to further verify the out-of-sample performance. It is worth noting that the SAA model can be viewed as a special case of the DRO model by restricting {\overline{p}}_s={\underset{\_}{p}}_s={\displaystyle \frac{1}{S}} for any s\in \mathcal{S} and {ϵ}_{i,t}={\sigma }_{i,j,t}=0 for any i,j\in \mathcal{R} and t\in \mathcal{T} in the ambiguity set \mathbb{F}, as shown in Proposition 2. Moreover, the RO models include the formulations with box uncertainty sets and budget uncertainty sets. Notice that the supervision strategy h_{i,t} is a binary random variable, and in reality, each region could be supervised in any period. Accordingly, the uncertainty set for supervision strategy h_{i,t} is formulated as {\mathcal{U}}_{\mathrm{h}}=\{\mathbf{h}|{\underset{\_}{h}}_{i,t}\le h_{i,t}\le {\overline{h}}_{i,t},\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T}\}, and we set {\underset{\_}{h}}_{i,t}=0 and {\overline{h}}_{i,t}=1. For the uncertain travel demand D_{i,j,t}, the box uncertainty set is shown in Eq. (12).

(12)\begin{array}{c}{\mathcal{U}}_{\mathrm{D}}^{\mathrm{b}\mathrm{o}\mathrm{x}}=\left\{\mathbf{D}|{\underset{\_}{D}}_{i,j,t}\le D_{i,j,t}\le {\overline{D}}_{i,j,t},\mathrm{\forall }i,j\in \mathcal{R},t\in \mathcal{T}\right\}.\end{array}

The box uncertainty set {\mathcal{U}}_{\mathrm{D}}^{\mathrm{b}\mathrm{o}\mathrm{x}} restricts the maximal and minimal value of the random demand D_{i,j,t}, which can be obtained from the historical data. However, the RO model with the box uncertainty set may lead to over-conservative results (Bertsimas and Sim, 2004). Moreover, we introduce the frequently used budget uncertainty set Eq. (13) as follows. Similar to Fu et al. (2022b), the budget uncertainty set (13) additionally includes the temporal correlations of O-D pairs to avoid over-conservative results, where {\mathrm{\Omega }}_t represents the maximal demand flow in each period t among all historical samples.

(13)\begin{array}{c}\begin{array}{c}{\mathcal{U}}_{\mathrm{D}}^{\mathrm{b}\mathrm{u}\mathrm{d}\mathrm{g}\mathrm{e}\mathrm{t}}=\{\mathbf{D}|{\underset{\_}{D}}_{i,j,t}\le D_{i,j,t}\le {\overline{D}}_{i,j,t},\\ \mathrm{\forall }i,j\in \mathcal{R},t\in \mathcal{T},\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}{D}_{i,j,t}\le {\mathrm{\Omega }}_t,\mathrm{\forall }t\in \mathcal{T}\}.\end{array}\end{array}

For the ambiguity set \mathbb{F} of the DRO model, {\mu }_{i,t} and {ϵ}_{i,t} are obtained according to the mean and mean absolute deviation of the generated supervision data. Similarly, {\theta }_{i,j,t} and {\sigma }_{i,j,t} are calculated according to the mean and mean absolute deviation of the historical demand data. In addition, the ranging parameters {\underset{\_}{p}}_s and {\overline{p}}_s that limit p_s are set to 0.9 and 1.1 respectively.

To conduct the simulation and verify the performance of the tested models, we randomly generated 2000 testing samples from the distribution estimated from historical data. Note that only the first-stage solutions are employed, and we re-optimize the second-stage decisions for each testing sample in the simulation procedure.

5.2 Out-of-sample model performance

We first assess the out-of-sample performance of the DRO model by comparing it with the benchmark models. We use several metrics, i.e., average value, 99th quantile, 95% CVaR and 99% CVaR, to assess the performances. The detailed definition of these metrics is shown in Table 2. Specifically, the average value reflects the out-of-sample average-case performance, and the other three metrics reflect the robust performance of the tested models.

Fig. 3 illustrates the out-of-sample performance of the DM, SAA, and DRO models. The deterministic model (DM) serves as the baseline, and the y-axis labeled “improved value” represents the difference in objective values between the DRO, SAA models and DM, specifically, the profit achieved by the DRO (SAA) solution minus that of the DM. A higher positive value indicates that the tested model yields better performance under the corresponding evaluation metric.

Fig. 3(a) shows the average value of the tested models with various bike budgets W. Compared with the DM, both the DRO and SAA models show much better performances. For the three robustness metrics—99th quantile, 99% CVaR, and 95% CVaR, the DRO and SAA models also significantly surpass the DM, which exhibit stronger out-of-sample robustness. We also notice that the DRO model outperforms the SAA model in terms of both average value and robustness metrics, highlighting the importance of accounting for the ambiguity in travel demand and supervision strategies.

Moreover, we test the out-of-sample performances of the DRO model compared with the RO models with box and budget uncertainty sets and show the results in Fig. 4. The results of the RO model with box uncertainty set are selected as the benchmark. We notice that the RO model with budget uncertainty set outperforms the RO model with box uncertainty set, except for the case when W=1400, which verifies the price of robustness (Bertsimas and Sim, 2004). In addition, the DRO model has significantly better performance than the two RO models. The reason is that the DRO model additionally includes distributional information like moment and historical sample data, thus achieving better out-of-sample performance. Besides, notice that the supervision strategy h_{i,t} is binary uncertain variable, and the RO models can only consider the worst-case supervision strategy, i.e., each region will be supervised, and accordingly the results will be over-conservative. These comparisons indicate that the RO models are not suitable to handle such problems, while the DRO model could address this drawback by incorporating the scenario information.

Furthermore, we also explore the effect of the ranging parameter {\underset{\_}{p}}_s and {\overline{p}}_s in the ambiguity set. We introduce an auxiliary parameter \delta and let {\underset{\_}{p}}_s=\left(1-\delta \right){\displaystyle \frac{1}{S}}, {\overline{p}}_s=\left(1+\delta \right){\displaystyle \frac{1}{S}}. We fix the budget of bikes W=1700 and let \delta \in \{0,0.05,0.1,0.15,0.2,0.25\} and test the effect of various \delta on the performance of the DRO model, and show the results in Fig. 5. Specifically, \delta =0 represents that the DRO model is reduced to the SAA model with empirical distribution. We use the average value and 95% CVaR metrics to represent the performance. The results suggest that the DRO model achieves the best out-of-sample performance when \delta =0.1, which indicates that the DRO model outperforms the SAA model (\delta =0) by incorporating the moment distributional information and allowing the variation of probability of each scenario. Furthermore, as \delta increases beyond 0.1, both the average value and the 95% CVaR decline. This may be attributed to the fact that larger \delta values lead to overly conservative solutions, which result in worse performance compared with the DRO model with a smaller \delta.

In addition, Table 3 reports the computational times of the tested models under different values of W. As expected, the deterministic model (DM) and RO model with box uncertainty set exhibit the shortest computational time, while the SAA, RO model with budget uncertainty set, and the DRO models require more time due to their increased complexity. Nevertheless, the average computational time of the DRO model remains acceptable for a strategic-level optimization problem.

5.3 Effect analysis and managerial insights

In this part, we propose several important research questions and aim to provide managerial insights for both the operators of bike-sharing systems and the government supervisors by answering these questions. Specifically, we focus on the following five questions:

Firstly, does the supervision strategy help mitigate the oversupplying of bikes? To answer this question, we first focus on the effect of the supervision strategy on profits and bike allocation decisions. In Fig. 6(a), we show the average value and 95% CVaR of the profit. Compared with the case without supervision, the supervision of the government will reduce the profit obtained by the operator, because the supervision strategy will restrict the maximum bike supply in one region, and the oversupply of bikes will result in the penalty cost. In addition, we illustrate the top 5 regions with the largest number of initially allocated bikes for the cases with and without supervision in Fig. 6(b). We note that the supervision strategy can significantly mitigate the oversupply of bikes compared with the case without supervision, which may contribute to improving the congestion caused by crowded bikes.

To further explore the effect of the supervision strategy, we show the probability and degree of oversupplying that occurs during the whole operational period \mathcal{T} among all regions \mathcal{R} in the out-of-sample test in Fig. 7. Specifically, the probability of oversupplying is defined by the number of regions where oversupplying occurs during the whole operational period \mathcal{T} divided by the total number of regions times operational periods (R\times T) in all testing samples. The degree of oversupply is calculated by the total number of oversupplying bikes divided by the frequency at which oversupplying occurs in all simulation samples. A lower value of oversupplying probability and degree means better operational order. In Fig. 7(a), we show the probability of oversupplying. Compared with the case without supervision, the supervision strategy can reduce the probability of oversupplying by 6%. Fig. 7(b) shows the degree of oversupplying, and we note that the supervision strategy also contributes to reducing the degree of oversupplying. Thus, the results suggest that the supervision strategy is efficient in reducing both the probability and degree of oversupplying of bikes and helps to build the operational order in the bike-sharing market.

Secondly, how does the maximum inventory levelUaffect the operation? To address this problem, we test the effect of the maximum inventory level U stipulated by the government on the performance and show the results in Fig. 8. Specifically, we vary the scale control parameter U_{\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{e}\mathrm{n}\mathrm{t}} from 0.8 to 1.2 with a step of 0.1, and set U_{i,t}=U_{\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{e}\mathrm{n}\mathrm{t}}\times U_{i,t}^{\mathrm{b}\mathrm{a}\mathrm{s}\mathrm{i}\mathrm{c}}. A larger U_{\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{e}\mathrm{n}\mathrm{t}} indicates that the bike-sharing system can allocate more bikes in the operational regions. In Fig. 8(a), we first illustrate the average value and the 95% CVaR of the profit gained by the operator. The results show that with the increase of U_{\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{e}\mathrm{n}\mathrm{t}}, the profit gained by the operator also increases because a larger U_{\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{e}\mathrm{n}\mathrm{t}} means that the operator can allocate more bikes in a certain region to satisfy more demand without paying the penalty for oversupplying.

Moreover, we show the effect of U_{\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{e}\mathrm{n}\mathrm{t}} on the initial bike allocation. In Fig. 8(b), we show the top 5 regions with the largest number of initially allocated bikes with various U_{\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{e}\mathrm{n}\mathrm{t}}. We notice that compared with the cases of higher U_{\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{e}\mathrm{n}\mathrm{t}}, the bike-sharing system tends to allocate more bikes in one region to satisfy more demands when U_{\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{e}\mathrm{n}\mathrm{t}}=0.8 and 0.9, which may violate the maximum inventory restriction U_{i,t} and incur penalty cost. The reason for such results may be that the travel demand in certain regions is significantly higher than in other regions, which highlights the necessity of considering regional heterogeneous characteristics when setting the maximum bike inventory. While the bike allocations in other regions are significantly reduced due to the restriction of a lower U_{\mathrm{c}\mathrm{o}\mathrm{e}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}\mathrm{i}\mathrm{e}\mathrm{n}\mathrm{t}}.

Thus, we suggest that the government supervisors should consider the demand pattern and the heterogeneous characteristics of each region, such as the density of population and travel demand among different operational regions, when determining the maximum inventory level U_{i,t}. A lower U_{i,t} may result in more unsatisfied demands and a heavy operational burden for the operators, while a higher U_{i,t} may not achieve the purpose of restricting the bike allocation within a reasonable range to mitigate the traffic congestion caused by crowded bikes.

Thirdly, how does the malfunction probability of available bikes\alphaaffect the operation? To answer this question, we focus on the operational performance with various malfunction probabilities of bikes \alpha. We vary \alpha \in \{0.01,0.02,0.03,0.04,0.05\} and illustrate the effect of various \alpha in Fig. 9. Specifically, we show the effect of \alpha on the profit obtained by the operator in Fig. 9(a). Intuitively, both the average value and the 95% CVaR of the profit decrease when \alpha gets larger. Furthermore, in Fig. 9(b), we show the number of broken bikes transported to depots during the whole day, and the result can explain the trend shown in Fig. 9(a). A larger \alpha indicates that more malfunctions of bikes will occur. In this case, there are fewer available bikes combined with the increased penalty cost for malfunctioning bikes, which will result in a decline in profit. Thus, we suggest that the operators should focus on the daily maintenance of the bikes and avoid frequent malfunctions of the bikes to improve the operational performance and profit.

Fourthly, do the bike rebalancing operation and the construction of depots improve the operation performance? To figure out this question, we show the effect of bike rebalancing operation and depots on the operation performance of the bike-sharing system in Fig. 10. Specifically, we compare the profit obtained by the operator in the case with bike rebalancing and depot, and the case without bike rebalancing and depot. The profit under the case with bike rebalancing and depot is significantly higher compared to the case without bike rebalancing and depot, in terms of both average value and 95% CVaR, which verifies the necessity of the bike rebalancing operation and depot. The reason is that the operator can move malfunctioning bikes to the depots through bike rebalancing to avoid the penalty caused by malfunctioning bikes, and relocate the oversupplying bikes to other regions to avoid the penalty for oversupply.

Lastly, how does the demand elasticity affect the bike-sharing system? Similar to Jiang et al. (2020) and Fu et al. (2022a), we assume that the demands of users are related to the bike supply level W, which satisfies the following relationship, where b_1 and b_2 represent the maximal total demand and the demand elasticity, respectively.

(14)\begin{array}{c}{\sum }_{i\in \mathcal{R}}{\sum }_{j\in \mathcal{R}}{\sum }_{t\in \mathcal{T}}{D}_{ijt}=b_1\mathrm{e}\mathrm{x}\mathrm{p}\left(\frac{-b_2}{\sqrt{W}}\right).\end{array}

Specifically, we test the effect of demand elasticity by setting b_2=10 and 50, to represent low demand elasticity and high demand elasticity, respectively. Then we use the case W=1300 and the historical demand as the baseline, and calculate the maximal total demand b_1 under low demand elasticity and high demand elasticity cases. We show how the demand elasticity affects the average profit of the bike-sharing system in Fig. 11. The results illustrate that the demand elasticity has little effect on the average profit when the bike supply level W is low, for example, when W=1300 and W=1400. However, when the operator of the bike-sharing system further increases the bike supply W, we notice that the profit also increases in the high demand elasticity case than in the low demand elasticity case. A potential reason is that when the bike supply W is low, there are not sufficient available bikes to serve the additional demand. When W increases, both the available and demand increase, thus leading to higher operational profit. Accordingly, we suggest that the operator can increase the bike budget when the demand elasticity is high in bike-sharing systems, which could improve the operational profit.

In practice, the operators of bike-sharing systems can implement the following steps to optimize the strategic and operational decisions. I. Data collection. Firstly, before making the strategic and operational decisions, the operators should collect the historical travel demand for the bike-sharing system in the planning area. The operators can use their own historical trip records as the dataset. II. Parameters calibration. Some important parameters in the model should be calibrated according to the real-world setting or simulation results, such as the maximum inventory restriction U_{i,t} and malfunction rate of bikes \alpha. Calibrating these parameters could guarantee the effectiveness of the optimization model and improve the real-world performance. III. Model solution. After data collection and parameter calibration, the operators could formulate and solve the proposed models, and obtain the related strategic-level decisions (e.g., depot location and bike allocation) and operational-level decisions (e.g., daily fleet operations). IV. Decision implementation. The strategic-level decisions could be implemented because these decisions are long-term fixed decisions. For the daily operation, the operators could observe the daily demand and supervision strategy for a certain day, and then re-optimize the proposed model (like the simulation procedure in the numerical experiments) to derive and implement the operational-level decisions, which follow an adaptive strategy.

6 Conclusions

This study addresses the integrated optimization of depot location, inventory management and bike rebalancing problems in bike-sharing systems, considering government supervision. The supervision for oversupplying bikes and malfunctioning bikes is considered. Incorporating the uncertainty of demand and government supervision strategy, we develop a DRO model to address uncertainty. Specifically, we determine the depot location and the bike allocation among the regions and depots in the first stage. The dynamic bike rebalancing operation with an inventory transshipment framework is formulated in the second stage.

The performances of the DRO model are tested via real-world data from Mobike. First, we verified that the DRO model outperforms the benchmark models in terms of a series of metrics. We also conduct an effect analysis to assess the performance of the government supervision strategy. The results suggest that the supervision strategy helps mitigate the oversupplying of bikes, and the government supervisors should consider the demand pattern in the operational regions when making the maximum inventory level decision. In addition, we also reveal that the malfunction probability of available bikes significantly affects the operational performance and suggest that the operators should focus on the daily maintenance of bikes. Furthermore, we verify the positive impact of bike rebalancing operations and depots in the bike-sharing system. Considering the demand elasticity, the bike-sharing system could increase the bike supply to obtain more profit.

There are some limitations in the current study. For example, this paper only considers the operational strategy of one bike-sharing system under government supervision. In fact, in a competitive bike-sharing market, the government supervision strategy also affects other companies. The operational strategy of one bike-sharing system is also affected by the changed operational strategies of other competitive companies. Thus, future research can consider the impact of government supervision in a competitive bike-sharing market. Additionally, a customized solution algorithm can be developed to address the proposed reformulated DRO model with larger sizes.

7 Appendixes

Proof of Proposition 1.

The inner minimization problem {\mathrm{i}\mathrm{n}\mathrm{f}}_{\mathbb{P}\in \mathbb{F}}{\mathbb{E}}_{\mathbb{P}}\left[H\left(\mathit{\nu },{\mathbf{h}}_{\mathbf{s}},{\mathbf{D}}_{\mathbf{s}}\right)\right] in Model (10) is represented as follows, given \mathbb{F} and \mathit{\nu }:

(A.1a)\begin{array}{c}\mathrm{m}\mathrm{i}\mathrm{n}\sum _{s\in \mathcal{S}}{p}_{s}H\left(\mathit{\nu },{\mathbf{h}}_{\mathbf{s}},{\mathbf{D}}_{\mathbf{s}}\right),\end{array}

(A.1b)\begin{array}{c}\text{s.t.}\sum _{s\in \mathcal{S}}{p}_s=1,\end{array}

(A.1c)\begin{array}{c}-\sum _{s\in \mathcal{S}}{h}_{i,t}^s{p}_s\ge -{\mu }_{i,t}-{ϵ}_{i,t},\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},\end{array}

(A.1d)\begin{array}{c}\sum _{s\in \mathcal{S}}{h}_{i,t}^s{p}_s\ge {\mu }_{i,t}-{ϵ}_{i,t},\mathrm{\forall }i\in \mathcal{R},t\in \mathcal{T},\end{array}

(A.1e)\begin{array}{c}-\sum _{s\in \mathcal{S}}{D}_{i,j,t}^s{p}_s\ge -{\theta }_{i,j,t}-{\sigma }_{i,j,t},\mathrm{\forall }i,j\in \mathcal{R},t\in \mathcal{T},\end{array}

(A.1f)\begin{array}{c}\sum _{s\in \mathcal{S}}{D}_{i,j,t}^s{p}_s\ge {\theta }_{i,j,t}-{\sigma }_{i,j,t},\mathrm{\forall }i,j\in \mathcal{R},t\in \mathcal{T},\end{array}

(A.1g)\begin{array}{c}{p}_s\ge {\underset{\_}{p}}_s,\mathrm{\forall }s\in \mathcal{S},\end{array}

(A.1h)\begin{array}{c}-p_s\ge -{\overline{p}}_s,\mathrm{\forall }s\in \mathcal{S},\end{array}

(A.1i)\begin{array}{c}{p}_s\ge 0,\mathrm{\forall }s\in \mathcal{S}.\end{array}

We denote \lambda, {\mathit{\eta }}^{\mathbf{1}}, {\mathit{\eta }}^{\mathbf{2}}, {\mathit{\varphi }}^{\mathbf{1}}, {\mathit{\varphi }}^{\mathbf{2}}, {\mathit{\psi }}^{\mathbf{1}}, {\mathit{\psi }}^{\mathbf{2}} as the dual variables of Constraints (A.1b)–(A.1h). Based on duality theory, odel (A.1) can be reformulated to maximization Model (A.2).

(A.2a)\begin{array}{c}\begin{array}{rl}& \mathrm{m}\mathrm{a}\mathrm{x}\lambda -\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}\left({\mu }_{i,t}+{ϵ}_{i,t}\right){\eta }_{i,t}^1+\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}\left({\mu }_{i,t}-{ϵ}_{i,t}\right){\eta }_{i,t}^2\\ & -\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}\left({\theta }_{i,j,t}+{\sigma }_{i,j,t}\right){\varphi }_{i,j,t}^1\\ & +\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}\left({\theta }_{i,j,t}-{\sigma }_{i,j,t}\right){\varphi }_{i,j,t}^2+\sum _{s\in \mathcal{S}}\left({\underset{\_}{p}}_s{\psi }_s^1-{\overline{p}}_s{\psi }_s^2\right),\end{array}\end{array}

(A.2b)\begin{array}{c}\begin{array}{rl}& \text{s.t.}\lambda -\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{h}_{i,t}^s{\eta }_{i,t}^1+\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{h}_{i,t}^s{\eta }_{i,t}^2-\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{D}_{i,j,t}^s{\varphi }_{i,j,t}^1\\ & +\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{D}_{i,j,t}^s{\varphi }_{i,j,t}^2+{\psi }_s^1-{\psi }_s^2\le H\left(\mathit{\nu },{\mathbf{h}}^{\mathbf{s}},{\mathbf{D}}^{\mathbf{s}}\right),\mathrm{\forall }s\in \mathcal{S},\end{array}\end{array}

(A.2c)\begin{array}{c}{\eta }_{i,t}^1,{\eta }_{i,t}^2,{\varphi }_{i,j,t}^1,{\varphi }_{i,j,t}^2,{\psi }_s^1,{\psi }_s^2\ge 0,\mathrm{\forall }s,\in \mathcal{S},i,j\in \mathcal{R},t\in \mathcal{T}.\end{array}

Incorporating the maximization of first-stage decisions \mathit{\nu }, we can obtain the reformulation shown in Proposition 1.

Proof of Proposition 2.

Suppose that {\overline{p}}_s={\underset{\_}{p}}_s={\displaystyle \frac{1}{S}} for any s\in \mathcal{S} and {ϵ}_{i,t}={\sigma }_{i,j,t}=0 for any i,j\in \mathcal{R} and t\in \mathcal{T} in the ambiguity set \mathbb{F}, the reformulated model of the DRO model is shown as follows:

(A.3a)\begin{array}{c}\begin{array}{rl}& \mathrm{m}\mathrm{a}\mathrm{x}\lambda -\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{\mu }_{i,t}\left({\eta }_{i,t}^1-{\eta }_{i,t}^2\right)-\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{\theta }_{i,j,t}\left({\varphi }_{i,j,t}^1-{\varphi }_{i,j,t}^2\right)\\ & +\sum _{s\in \mathcal{S}}\frac{1}{S}\left({\psi }_s^1-{\psi }_s^2\right),\end{array}\end{array}

(A.3b)\begin{array}{c}\begin{array}{rl}& \text{s.t.}\lambda -\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{h}_{i,t}^s{\eta }_{i,t}^1+\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{h}_{i,t}^s{\eta }_{i,t}^2-\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{D}_{i,j,t}^s{\varphi }_{i,j,t}^1\\ & +\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{D}_{i,j,t}^s{\varphi }_{i,j,t}^2+{\psi }_s^1-{\psi }_s^2=H\left(\mathit{\nu },{\mathbf{h}}^{\mathbf{s}},{\mathbf{D}}^{\mathbf{s}}\right),\\ & \mathrm{\forall }s\in \mathcal{S},\end{array}\end{array}

(A.3c)\begin{array}{c}{\eta }_{i,t}^1,{\eta }_{i,t}^2,{\phi }_{i,j,t}^1,{\phi }_{i,j,t}^2,{\psi }_s^1,{\psi }_s^2\ge 0,\mathrm{\forall }s,\in \mathcal{S},i,j\in \mathcal{R},t\in \mathcal{T}.\end{array}

In Model (A.3), Constraint (A.3b) is an equation because p_s={\displaystyle \frac{1}{S}} and we omit p_s\ge 0 for any s\in \mathcal{S}. We introduce an auxiliary variable \mathrm{\Gamma } and reformulate Model (A.3) as follows:

(A.4a)\begin{array}{c}\mathrm{m}\mathrm{a}\mathrm{x}\mathrm{\Gamma },\end{array}

(A.4b)\begin{array}{c}\begin{array}{rl}& \text{s.t.}\mathrm{\Gamma }\le \lambda -\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{\mu }_{i,t}\left({\eta }_{i,t}^1-{\eta }_{i,t}^2\right)\\ & -\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{\theta }_{i,j,t}\left({\varphi }_{i,j,t}^1-{\varphi }_{i,j,t}^2\right)+\sum _{s\in \mathcal{S}}\frac{1}{S}\left({\psi }_s^1-{\psi }_s^2\right),\end{array}\end{array}

(A.4c)\begin{array}{c}\begin{array}{rl}& \lambda -\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{h}_{i,t}^s{\eta }_{i,t}^1+\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{h}_{i,t}^s{\eta }_{i,t}^2-\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{D}_{i,j,t}^s{\varphi }_{i,j,t}^1\\ & +\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{D}_{i,j,t}^s{\varphi }_{i,j,t}^2+{\psi }_s^1-{\psi }_s^2=H\left(\nu ,{\mathit{h}}^{\mathit{s}},{\mathit{D}}^{\mathit{s}}\right),\mathrm{\forall }s\in \mathcal{S},\end{array}\end{array}

(A.4d)\begin{array}{c}{\eta }_{i,t}^1,{\eta }_{i,t}^2,{\phi }_{i,j,t}^1,{\phi }_{i,j,t}^2,{\psi }_s^1,{\psi }_s^2\ge 0,\mathrm{\forall }s\in \mathcal{S},i,j\in \mathcal{R},t\in \mathcal{T}.\end{array}

For Constraint (A.4c), let both the left-hand-side and the right-hand-side formulation times {\displaystyle \frac{1}{S}}, then we note that the right-hand-side of Constraint (A.4b) is exactly the summation of the left-hand-side of Constraint (A.4c) for all s\in \mathcal{S} because {\mu }_{i,t}={\displaystyle \frac{\sum _{s\in \mathcal{S}}{h}_{i,t}^s}{S}} and {\theta }_{i,j,t}={\displaystyle \frac{\sum _{s\in \mathcal{S}}{D}_{i,j,t}^s}{S}}. Thus, we can rewrite Model (A.4) into following formulation:

(A.5a)\begin{array}{c}\mathrm{m}\mathrm{a}\mathrm{x}\mathrm{\Gamma },\end{array}

(A.5b)\begin{array}{c}\text{s.t.}\mathrm{\Gamma }\le \sum _{s\in \mathcal{S}}\frac{1}{S}H\left(\mathit{\nu },{\mathbf{h}}^{\mathbf{s}},{\mathbf{D}}^{\mathbf{s}}\right).\end{array}

Removing the auxiliary variable \mathrm{\Gamma } and incorporating the constraints for first-stage decisions \nu, we note that Model (A.5) is equivalent to the stochastic optimization model proposed in Section 4.1 and reach the conclusion in Proposition 2.

Formulation of DRO model

The detailed formulation of the reformulated equivalent deterministic model of the DRO model is shown as follows:

(B.1a)\begin{array}{c}\begin{array}{rl}& \mathrm{m}\mathrm{a}\mathrm{x}\lambda -\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}\left({\mu }_{i,t}+{ϵ}_{i,t}\right){\eta }_{i,t}^1+\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}\left({\mu }_{i,t}-{ϵ}_{i,t}\right){\eta }_{i,t}^2\\ & -\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}\left({\theta }_{i,j,t}+{\sigma }_{i,j,t}\right){\varphi }_{i,j,t}^1+\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}\left({\theta }_{i,j,t}-{\sigma }_{i,j,t}\right){\varphi }_{i,j,t}^2\\ & +\sum _{s\in \mathcal{S}}\left({\underset{\_}{p}}_s{\psi }_s^1-{\overline{p}}_s{\psi }_s^2\right),\end{array}\end{array}

(B.1b)\begin{array}{c}\begin{array}{rl}& \text{s.t.}\lambda -\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{h}_{i,t}^s{\eta }_{i,t}^1+\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{h}_{i,t}^s{\eta }_{i,t}^2-\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{D}_{i,j,t}^s{\varphi }_{i,j,t}^1\\ & +\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{D}_{i,j,t}^s{\varphi }_{i,j,t}^2+{\psi }_s^1-{\psi }_s^2\le \sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}r{d}_{i,j,t}^s\\ & -\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{c}_1{u}_{i,j,t}^s-\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{c}_2{v}_{i,t}^s-\sum _{i\in \mathcal{R}}\sum _{t\in \mathcal{T}}{c}_3{n}_{i,t}^{\mathrm{b},s}\\ & -\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\sum _{t\in \mathcal{T}}{c}_{4,i,j}\left(k_{i,j,t}^s+y_{i,j,t}^s+z_{i,j,t}^s\right),\mathrm{\forall }s\in \mathcal{S},\end{array}\end{array}

(B.1c)\begin{array}{c}\sum _{i\in \mathcal{R}}{x}_i=P,\end{array}

(B.1d)\begin{array}{c}\sum _{i\in \mathcal{R}}N{n}_i^a+\sum _{i\in \mathcal{R}}N{q}_i^a=W,\end{array}

(B.1e)\begin{array}{c}N{q}_i^a\le M{x}_i,\mathrm{\forall }i\in \mathcal{R},\end{array}

(B.1f)\begin{array}{c}N{n}_i^a=n_{i,1}^{a,s},\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},\end{array}

(B.1g)\begin{array}{c}N{q}_i^a=q_{i,1}^{a,s},\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},\end{array}

(B.1h)\begin{array}{c}{n}_{i,1}^{b,s}=0,\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},\end{array}

(B.1i)\begin{array}{c}{q}_{i,1}^{b,s}=0,\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},\end{array}

(B.1j)\begin{array}{c}{q}_{i,t}^{a,s}\le Q^a{x}_i,\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},t\in \mathcal{T},\end{array}

(B.1k)\begin{array}{c}{q}_{i,t}^{b,s}\le Q^b{x}_i,\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},t\in \mathcal{T},\end{array}

(B.1l)\begin{array}{c}\sum _{j\in \mathcal{R}}\left(d_{i,j,t}^s+k_{i,j,t}^s\right)\le n_{i,t}^{a,s},\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},t\in \mathcal{T},\end{array}

(B.1m)\begin{array}{c}\sum _{j\in \mathcal{R}}{z}_{i,j,t}^s\le n_{i,t}^{b,s},\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},t\in \mathcal{T},\end{array}

(B.1n)\begin{array}{c}\sum _{j\in \mathcal{R}}{y}_{i,j,t}^s\le q_{i,t}^{a,s},\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},t\in \mathcal{T},\end{array}

(B.1o)\begin{array}{c}{y}_{i,j,t}^s\le M{x}_i,\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},\end{array}

(B.1p)\begin{array}{c}{\mathit{\text{z}}}_{j,i,t}^s\le M{x}_i,\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},\end{array}

(B.1q)\begin{array}{c}{d}_{i,j,t}^s\le D_{i,j,t}^s,\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},\end{array}

(B.1r)\begin{array}{c}{u}_{i,j,t}^s\ge D_{i,j,t}^s-d_{i,j,t}^s,\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},\end{array}

(B.1s)\begin{array}{c}\sum _{i\in \mathcal{R}}\sum _{j\in \mathcal{R}}\left(k_{i,j,t}^s+y_{i,j,t}^s+{\mathit{\text{z}}}_{i,j,t}^s\right)\le C,\mathrm{\forall }s\in \mathcal{S},t\in \mathcal{T},\end{array}

(B.1t)\begin{array}{c}{e}_{i,t}^s\ge n_{i,t}^{a,s}+n_{i,t}^{b,s}-U_{i,t},\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},t\in \mathcal{T},\end{array}

(B.1u)\begin{array}{c}{v}_{i,t}^s\ge h_{i,t}^s{e}_{i,t}^s,\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},t\in \mathcal{T},\end{array}

(B.1v)\begin{array}{c}{q}_{i,t}^{a,s}=q_{i,t-1}^{a,s}-\sum _{j\in \mathcal{R}}{y}_{i,j,t-1}^s,\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},t\in \mathcal{T},t\ge 1,\end{array}

(B.1w)\begin{array}{c}{q}_{i,t}^{b,s}=q_{i,t-1}^{b,s}+\sum _{j\in \mathcal{R}}{\mathit{\text{z}}}_{j,i,t-1}^s,\mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},t\in \mathcal{T},t\ge 1,\end{array}

(B.1x)\begin{array}{c}\begin{array}{rl}{n}_{i,t}^{a,s}=& n_{i,t-1}^{a,s}-\sum _{j\in \mathcal{R}}\left(d_{i,j,t-1}^s+k_{i,j,t-1}^s\right)\\ & +\sum _{j\in \mathcal{R}}\left(\left(1-\alpha \right)d_{j,i,t-1}^s+k_{j,i,t-1}^s+y_{j,i,t-1}^s\right),\\ & \mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},t\in \mathcal{T},t\ge 1,\end{array}\end{array}

(B.1y)\begin{array}{c}\begin{array}{r}{n}_{i,t}^{b,s}=n_{i,t-1}^{b,s}-\sum _{j\in \mathcal{R}}{z}_{i,j,t-1}^s+\alpha \sum _{j\in \mathcal{R}}{d}_{j,j,t-1}^s,\\ \mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},t\in \mathcal{T},t\ge 1,\end{array}\end{array}

(B.1z)\begin{array}{c}{x}_i\in \{0,1\},\mathrm{\forall }i\in \mathcal{R},\end{array}

(B.1aa)\begin{array}{c}N{n}_i^a,N{q}_i^a\in {\mathbb{Z}}^{+},\mathrm{\forall }i\in \mathcal{R},\end{array}

(B.1ab)\begin{array}{c}\begin{array}{c}{q}_{i,t}^{a,s},q_{i,t}^{b,s},n_{i,t}^{a,s},n_{i,t}^{b.s},e_{i,t}^s,v_{i,t}^s,k_{i,j,t}^s,d_{i,j,t}^s,{\mathit{\text{z}}}_{i,j,t}^s,y_{i,j,t}^s,u_{i,j,t}^s\in {\mathbb{Z}}^{+},\\ \mathrm{\forall }s\in \mathcal{S},i\in \mathcal{R},j\in \mathcal{R},t\in \mathcal{T},\end{array}\end{array}

(B.1ac)\begin{array}{c}{\eta }_{i,t}^1,{\eta }_{i,t}^2,{\phi }_{i,j,t}^1,{\phi }_{i,j,t}^2,{\psi }_s^1,{\psi }_s^2\ge 0,\mathrm{\forall }s\in \mathcal{S},i,j\in \mathcal{R},t\in \mathcal{T}.\end{array}

References

[1] Ben-Tal A, Den Hertog D, De Waegenaere A, Melenberg B, Rennen G, (2013). Robust solutions of optimization problems affected by uncertain probabilities. Manage Sci, 59: 341–357

[2] Ben-Tal A, Nemirovski A, (1999). Robust solutions of uncertain linear programs. Oper Res Lett, 25: 1–13

[3] Bertsimas D, Sim M, (2004). The price of robustness. Oper Res, 52: 35–53

[4] Birge J RLouveaux F (2011). Introduction to stochastic programming. Springer Science & Business Media

[5] Bulhões T, Subramanian A, Erdogan G, Laporte G, (2018). The static bike relocation problem with multiple vehicles and visits. Eur J Oper Res, 264: 508–523

[6] Caggiani L, Colovic A, Ottomanelli M, (2020). An equality-based model for bike-sharing stations location in bicycle-public transport multimodal mobility. Transp Res Part A Policy Pract, 140: 251–265

[7] Chang X, Wu J, Sun H, Yan X, (2023). A smart predict-then-optimize method for dynamic green bike relocation in the free-floating system. Transp Res, Part C Emerg Technol, 153: 104220

[8] Chen Q, Fu C, Zhu N, Ma S, He Q, (2023). A target-based optimization model for bike-sharing systems: From the perspective of service efficiency and equity. Transp Res, Part B: Methodol, 167: 235–260

[9] Datner S, Raviv T, Tzur M, Chemla D, (2019). Setting inventory levels in a bike sharing network. Transport Sci, 53: 62–76

[10] Erdogan G, Battarra M, Wolfler Calvo R, (2015). An exact algorithm for the static rebalancing problem arising in cicycle sharing systems. Eur J Oper Res, 245: 667–679

[11] Mohajerin Esfahani P, Kuhn D, (2018). Data-driven distributionally robust optimization using the wasserstein metric: Performance guarantees and tractable reformulations. Math Program, 171: 115–166

[12] Fishman E, (2016). Bikeshare: A review of recent literature. Transport Reviews, 36: 92–113

[13] Frade I, Ribeiro A, (2015). Bike-sharing stations: A maximal covering location approach. Transp Res Part A Policy Pract, 82: 216–227

[14] Fu C, Ma S, Zhu N, He Q, Yang H, (2022a). Bike-sharing inventory management for market expansion. Transp Res, Part B: Methodol, 162: 28–54

[15] Fu C, Zhu N, Ma S, Liu R, (2022b). A two-stage robust approach to integrated station location and rebalancing vehicle service design in bike-sharing systems. Eur J Oper Res, 298: 915–938

[16] Guo Y, Li J, Xiao L, Allaoui H, Choudhary A, Zhang L, (2024). Efficient inventory routing for Bike-Sharing Systems: A combinatorial reinforcement learning framework. Transp Res, Part E Logist Trans Rev, 182: 103415

[17] He L, Hu Z, Zhang M, (2020). Robust repositioning for vehicle sharing. Manuf Serv Oper Manag, 22: 241–256

[18] Jiang Z, Lei C, Ouyang Y, (2020). Optimal investment and management of shared bikes in a competitive market. Transp Res, Part B: Methodol, 135: 143–155

[19] Jin Y, Ruiz C, Liao H, (2022). A simulation framework for optimizing bike rebalancing and maintenance in large-scale bike-sharing systems. Simul Model Pract Theory, 115: 102422

[20] Kaspi M, Raviv T, a Tzur M, (2016). Detection of unusable bicycles in bike-sharing systems. Omega, 65: 10–16

[21] Li L, Li G, Liang S, (2021). Does government supervision suppress free-floating bike sharing development? Evidence from Mobike in China. Inf Technol Dev, 27: 802–826

[22] Li X, Wang X, Feng Z, (2024). Dynamic repositioning in bike-sharing systems with uncertain demand: An improved rolling horizon framework. Omega, 126: 103047

[23] Lin J, Yang T, (2011). Strategic design of public bicycle sharing systems with service level constraints. Transp Res, Part E Logist Trans Rev, 47: 284–294

[24] Lin J, Yang T, Chang Y, (2013). A hub location inventory model for bicycle sharing system design: Formulation and solution. Comput Ind Eng, 65: 77–86

[25] Liu S, Shen Z M, Ji X, (2022). Urban bike lane planning with bike trajectories: Models, algorithms, and a real-world case study. Manuf Serv Oper Manag, 24: 2500–2515

[26] Long D Z, Sim M, Zhou M, (2023). Robust satisficing. Oper Res, 71: 61–82

[27] Lu C, (2016). Robust multi-period fleet allocation models for bike-sharing systems. Netw Spat Econ, 16: 61–82

[28] Lv C, Liu Q, Zhang C, Ren Y, Zhou H, (2024). A feature correlation reinforce clustering and evolutionary algorithm for the green bike-sharing reposition problem. Comput Oper Res, 166: 106627

[29] Lv C, Zhang C, Lian K, Ren Y, Meng L, (2022). A two-echelon fuzzy clustering based heuristic for large-scale bike sharing repositioning problem. Transp Res, Part B: Methodol, 160: 54–75

[30] Maggioni F, Cagnolari M, Bertazzi L, Wallace SW, (2019). Stochastic optimization models for a bike-sharing problem with transshipment. Eur J Oper Res, 276: 272–283

[31] Ministry of Transport in Singapore (2022). Written reply to parliamentary question on measures taken to ensure non-dumping or illegal parking of shared bicycles from granting of new licence for bike sharing. Available at the website of mot.gov.sg

[32] Nikiforiadis A, Aifadopoulou G, Salanova Grau JM, Boufidis N, (2021). Determining the optimal locations for bike-sharing stations: methodological approach and application in the city of Thessaloniki, Greece. Transp Res Procedia, 52: 557–564

[33] Pu D, Xie F, Yuan G, (2020). Active supervision strategies of online ride-hailing based on the tripartite evolutionary game model. IEEE Access, 8: 149052–149064

[34] Qin M, Wang J, Chen WM, (2023). Reducing CO₂ emissions from the rebalancing operation of the bike-sharing system in beijing. Frontiers of Engineering Management, 10(2): 262–284

[35] Raviv T, Kolka O, (2013). Optimal inventory management of a bike-sharing station. IIE Trans, 45: 1077–1093

[36] Raviv T, Tzur M, Forma IA, (2013). Static repositioning in a bike-sharing system: models and solution approaches. EURO J Transp Logist, 2: 187–229

[37] Schuijbroek J, Hampshire R C, Van Hoeve W J, (2017). Inventory rebalancing and vehicle routing in bike sharing systems. Eur J Oper Res, 257: 992–1004

[38] Shenzhen Municipal Transportation Bureau (2017). Policy interpretation of the “implementation plan for the regulation and rectification of internet rental bicycles in shenzhen”. Available at the website of jtys.sz.gov.cn

[39] Shu J, Chou MC, Liu Q, Teo C, Wang I, (2013). Models for effective deployment and redistribution of bicycles within public bicycle-sharing systems. Oper Res, 61: 1346–1359

[40] Shui CS, Szeto W, (2018). Dynamic green bike repositioning problem—A hybrid rolling horizon artificial bee colony algorithm approach. Transp Res Part D Transp Environ, 60: 119–136

[41] Shui CS, Szeto WY, (2020). A review of bicycle-sharing service planning problems. Transp Res, Part C Emerg Technol, 117: 102648

[42] Silva M C M, Aloise D, Jena S D, (2024). Data-driven prioritization strategies for inventory rebalancing in bike-sharing systems. Omega, 129: 103141

[43] Smith J E, Winkler R L, (2006). The optimizer’s curse: Skepticism and postdecision surprise in decision analysis. Manage Sci, 52: 311–322

[44] Song J, Li B, Szeto W, Zhan X, (2024). A station location design problem in a bike-sharing system with both conventional and electric shared bikes considering bike users’ roaming delay costs. Transp Res, Part E Logist Trans Rev, 181: 103350

[45] Wang Y, Jin H, Zheng S, Shang W, Wang K, (2023). Bike-sharing duopoly competition under government regulation. Appl Energy, 343: 121121

[46] Wiesemann W, Kuhn D, Sim M, (2014). Distributionally robust convex optimization. Oper Res, 62: 1358–1376

[47] Wu L, Gu W, Fan W, Cassidy M J, (2020). Optimal design of transit networks fed by shared bikes. Transp Res, Part B: Methodol, 131: 63–83

[48] Xiamen Municipal Law Enforcement Bureau of Urban Management (2025). Public consultation on the “2025 implementation plan for the assessment of internet rental bicycle operating companies in xiamen (draft for public comment)”. Available at the website of cg.xm.gov.cn

[49] Xiao G, Wang Z, (2020). Empirical study on bikesharing brand selection in china in the post-sharing era. Sustainability (Basel), 12: 3125

[50] Yan S, Lin J, Chen Y, Xie F, (2017). Rental bike location and allocation under stochastic demands. Comput Ind Eng, 107: 1–11

[51] Yin Z, Hardaway K, Feng Y, Kou Z, Cai H, (2023). Understanding the demand predictability of bike share systems: A station-level analysis. Frontiers of Engineering Management, 10: 551–565

[52] Zhang D, Yu C, Desai J, Lau H, Srivathsan S, (2017). A time-space network flow approach to dynamic repositioning in bicycle sharing systems. Transp Res, Part B: Methodol, 103: 188–207

[53] Zhao SHou WWang ZHe Q (2022). Bikes detention game in free-float bike sharing systems. In: INFORMS International Conference on Service Science: 251–26

[54] Zhong Y, Yang T, Cao B, Cheng T, (2022). On-demand ride-hailing platforms in competition with the taxi industry: Pricing strategies and government supervision. Int J Prod Econ, 243: 108301

RIGHTS & PERMISSIONS

Higher Education Press