The Integration of Multimodal Networks: The Generalized Modal Split and Collaborative Optimization of Transportation Hubs

,


Introduction
Te rapid development of technologies has led to the emergence of various new transport modes and mobility services. Te multimodal transport network, which consisted of these new forms of urban mobility and traditional transport, plays an overriding role in urban development [1,2]. However, the spatial in-coordination between the demand and multimodal capacity leads to some negative efects, such as trafc congestion, air pollution, and social equality problems, that makes the system unsustainable [3]. Concerning these problems, cities have realized that the key issue to solve the problems is to encourage a modal shift from car to public transportation (PT) [4].
Consequently, integrated transport has received a lot of attention from policymakers and researchers in recent 30 years because of its potential to foster sustainable travel behaviours [5,6].
As most cities have an existing multimodal network, it can be time and cost consuming for governments to create the changes in the current network [7]. From the infrastructure perspective, there are two ways to promote the attractiveness of PT. One is to provide more direct routes to reduce transfers [8]. Te other is to connect the routes with one another to create access a wider range of destinations [9]. Extended to the entire network, for cities with mature multimodal networks, it would be cost-efective to integrate all the travel modes through hubs rather than building more PT routes [10]. Urban transportation hubs, such as parking lots, metro stations, bus stops, and bicycle stations, are the places where transfers occur. Tey play the key role in integrating the multimodal network [11,12]). Tus, this paper focuses particularly on integrating the physical multimodal transport network via the optimization of transportation hubs to achieve a sustainable modal shift.
In the multimodal network, private car and conventional PT services are often considered as the major travel modes, complemented by shared, on-demand services and other travel modes [13]. Tey all have their own optimal ranges and advantages. Te integration aims to provide the users with all the available travel options to develop a customizable door-to-door combined trip that could increase the attractiveness of PT [14]. Recently, mobility as a service (MaaS) is considered as an evolutionary continuation in terms of integrated transport, and it is often described as a personalised user-oriented platform that unifes the creation, purchase, and delivery of door-to-door passenger transport solutions integrating traditional and innovative modes of transport [2,15,16]. Similarly, the success of MaaS depends on a shift in behaviour away from reliance upon private car ownership [17].
However, the efects of the combined mobility services are often questioned by the recent research works and applications [18], that MaaS should be complementary rather than in substitution to a private car. Even if in principle it may reduce the car usage of MaaS users, it may imply bigger use of the private car by non-MaaS users because of the crowding of transport services by MaaS users [19]. Tis reminds us that the sustainable mobility not only relies on a modal shift from car to PT but also on a balance between the demand and multimodal capacity. Since the efect of the integration did not happen immediately but gradually over the years, comparing with the short-term modal shift stimulated by new policies or reforms, the balanced use between all travel modes can keep it remain motivated [20,21]. Tus, the aim of the integration via hubs in this paper is extended to achieve a reasonable modal shift to promote a balanced use of all travel modes.
From a user's perspective, they have many mode choices to fnish their trips, including single modes and combined modes. A combined trip is usually less competitive than a single-mode trip, considering the transfer penalty. A user will only choose a combined travel path when the congestion cost of a single-mode path is higher than the transfer cost of the combined trip [22]. Tus, integrating a multimodal network by optimizing the hubs is a process that provides users with competitive combined travel choices by reducing transfer penalties. Te mixed choices among all modes of all users can realize a total modal shift from congested subnetwork to underutilized subnetwork while cutting down the usage of the private car.
In summary, to achieve a modal shift from congested subnetworks to underutilized subnetworks and promote a balanced use of the entire network, this study conducted the integration of the multimodal network by collaboratively optimizing the hubs of diferent subnetworks. Tis problem is referred to as the integrated multimodal network design problem (IMNDP) in this paper, focusing on the two key attributes of IMNDP, the combined modes and optimization of hubs. A new bilevel mathematical formulation is proposed based on the network design problem in this paper. To handle with the numerous modes options the network can provide, a generalized modal/split trafc assignment problem (GMS/TAP) is proposed in the lower-level problem. Meanwhile, for integrating the multimodal network through a collaborative optimization of transportation hubs, the methodologies for the assessment and optimal design of a hub integration scheme are provided in the upper-level problem.
Te remainder of this paper is organized as follows. Te literature review and research gaps are introduced in Section 2. Section 3 introduces the formulation of the lower-level GMS/TAP and upper-level problem along with their solution methods. Section 4 verifes the solution properties of the lower-level and upper-level problems along with the application in a real-size network. Te conclusions are fnally presented in Section 5.

Literature Review
Most of the previous researches on the integrated transportation focus on the policy side. UK Government addressed integration as "integration within and between diferent types of transport so that each contributes its full potential and people can move easily between them" [23]. Hull [24] defned eight levels for an integrated transport policy: (8) the integration of policy measures, (7) the integration of policy sectors, (6) institutional and administrative integration, (5) the integration of environmental issues, (4) integration with social objectives, (3) integration with market needs, (2) modal integration, and (1) physical and operational integration. Te diferent elements involved, especially the policy and institutional integration from approximately level 4 to 8, are very difcult to achieve, that makes integrated transport to be described as an inefective practice [25].
Among the 8 integration levels, the integration levels 1-2 are achievable [26]. Chowdhury et al. [7] stated that transfer is the key component of an integrated system, and proposed fve main attributes of an integrated system: (1) network integration, (2) fare integration, (3) information integration, (4) the physical integration of stations, and (5) coordinated schedules. Chowdhury et al. [27] showed that physical integration is more important than information integration for current users. Ibrahim [11] proposed that the network integration is a fully multimodal integrated system that optimizes the available resources and therefore reduces the wasteful duplication of services while still maintaining an adequate coverage area. Studies have also shown that network integration can increase the usage of PT services. Matas [28] investigated the signifcant increase in PT use (>40%) in Madrid. Te study discussed that integrated fare system and network integration had the most impact on ridership. Hidalgo and King [29] discussed the impact of the gradual implementation of an integrated PT system, which requires a combination of political leadership, technical capacity, and funding.
Following the prior eforts of the integrated transportation, with the development of mobile Internet and smartphone apps, MaaS is considered as an emerging opportunity to integrate the conventional PT, new on-demand, and shared mobility services. Sochor et al. [30] put forward three major attributes of the descriptions: (1) Ofering a service with customer/user/traveller/ consumer transport needs as the main focus (2) Ofering (multimodal) mobility rather than transport (3) Ofering integration of transport services, information, payment, and ticketing Te way to ofer multimodal mobility services via a mobile app, not the transport modes itself, makes the combined trip more acceptable for the users. Te combination of transport modes with diferent strengths and weaknesses is believed to be more efcient than one single mode [18]. However, even if MaaS has received a lot of attention, MaaS applications are still not a fully fedged reality. Many problems prevent its implementation and make the gap between theory and practice still high. Polydoropoulou et al. [2] elaborated on possible MaaS barriers and enablers related to infrastructures (integrated payments and open application programming interfaces), hard institutions (regulations, e.g., related to fares and security concerns), soft institutions (acceptability and trust between operators), and capabilities (need for transport investments). Some authors also highlight that bundling of mobility services is typically seen as the next step after the full integration of operations, information, and payments [6] since the infrastructure is the foundation of MaaS, providing the physical means to get from the origin to destination. Te network integration via the optimization of hubs in this paper can be viewed as a hardware support for the MaaS applications.
Quantitative studies to facilitate the integration scheme on modal shift are relatively sparse. Practical feld tests are regularly used in the later phases of the development process, but in the earlier phases, one is generally dependent on less powerful methods, such as checklists and expert judgements [31]. Solecka [32] proposed a common variant ranking method, Electre III, to assist policy makers in deciding the most suitable version of an integrated PT system for their city. Kash and Hidalgo [33] developed a framework to include consultation with the users and incorporate their needs into a regional PT plan. Data were collected by undertaking semistructured interviews with transportation experts. Considering the rigorous amounts of data and analysis which consumes time and resources, El-Adaway et al. [34] used the social network analysis method as a complementary tool for improved transportation planning. Other studies focused on the simulation method, Santi et al. [35] proposed a method of shareability networks, which formulates the spatiotemporal sharing problem using a graph-theoretic framework. Simulations on a data set of New York City taxi trips indicate that the cumulative trip length can be reduced by 40% through ride sharing. Shen et al. [36] designed an integrated autonomous vehicles (AV)-PT system by repurposing low-demand bus routes using the shared AV service while retaining the high-demand bus routes. Trough the agent-based simulations, the authors calculate the best number of AV vehicles to replace bus routes.
Since the plan of the network infrastructure afects the long-term passenger distributions. Te simulation methods concentrate on the dynamic and short-term evolution of the user's behaviour, which is good at predicting the usage of ondemand or real-time trafcs. However, the backbone of the multimodal network is metro, scheduled bus and car, integrating the hubs mainly serves these subnetworks. Tus, the multimodal network design problem (MNDP) which uses static trafc assignment problem is more suitable for capturing the fnal equilibrium fow among the entire network.
Te MNDP is often formulated as a bilevel mathematical program with equilibrium constraints [37][38][39], at the upperlevel, where the optimization variables focus on [40] (1) route design [41], (2) frequency setting [42,43], (3) timetabling [44], (4) vehicle scheduling [45], and (5) crew scheduling. MNDP is also classifed according to the aimed design period into strategic, tactical, and operational stages. Route design constitutes the strategic stage, and the frequency setting is undertaken at the tactical stage, and the remaining activities are handled at the operational stage. Te lower-level problem is usually a trafc assignment problem (TAP). Te indexes generated by the TAP will be fed back to the upper-level problem to assess the optimization schemes.
As a widely used tool to predict the fow distribution of a network, Shef [46] described the TAP of capturing the mode choices of users as a modal split/trafc assignment problem (MS/TAP), in which the fow of each mode was equivalent to postulating a travel demand function, where the number of trips was proportional to a negative exponential function of the travel cost and linked to route choice models based on the Wardropian user equilibrium [47]. Te single-mode model was frst successfully formulated and solved by Evans [48]. Ten, Leblanc [49]; Safwat and Magnanti [50]; Boyce and Zhang [51]; and many others studied various [52][53][54] aspects of the problem, including its formulations, algorithms, cost functions, and calibration methods. Among these studies, mode choice was often made between two or three single modes (car and transit) and combined with the route choice to describe the users' spatial distribution in the network [55,56].
As the single-mode models were being perfected, the research attention gradually moved to combined modes [57][58][59]. Many studies about specifc combined modes, including P&R (park and ride) [60], remote P&R [61], and transit and e-hailing [62], were proposed. Te mode-choice structure then changes from a single level to multiple levels and is usually formulated with a nested logit (NL) or crossnested logit (CNL) model [63]).
Te upper-level part of the MNDP applies a process similar to that from a single network [64][65][66] to multiple networks. Only a few researchers have simultaneously tackled over two subnetwork facilities or services. Chen and Nie [62] developed two diferent relative spatial position designs for an integrated e-hailing/fxed-route transit system. Song et al. [67] proposed an integrated planning framework to locate park and ride (P&R) facilities and simultaneously optimize their capacities and transit service frequencies. Tese studies had the goal of improving the attractiveness of one combined trip by simultaneously optimizing the services of the two modes involved.
As discussed above, the IMNDP of this paper focuses on integrating the physical multimodal network via the optimization of transportation hubs to achieve a balanced fow distribution within the entire network. Tis brings three problems in adopting it into the MNDP model structure: (1) IMNDP involves the overall optimization of hubs in diferent subnetworks, where the change in one subnetwork's hubs will infuence the fow distribution among other subnetworks, so the mode choice in the MS/TAP should consider all the modes that the multimodal network can provide. In the traditional MS/TAP, the mode choice only covers the aimed design modes. However, with the development of multimodal networks, more basic travel modes are introduced in while the increase in travel distance also brings in an increase of transfers, which leads to the combined modes the multimodal network can provide are almost incalculable.
(2) Unlike the MNDP, which focuses on increasing one subnetwork's performance, the IMNDP concentrates on the interaction between diferent modes. Te objectives of the upper-level problem should refect the modal shift and balanced use of the entire network. (3) Considering the collaborative optimization of different subnetwork hubs, a large search space requires highly efcient algorithms to solve the bilevel problem.
For the frst problem, if we use the mode structure of the traditional MS/TAP, the large number of modes will make the mathematical formulation extremely complicated, which can be found in Liu et al. [61]. In addition, for a path-based algorithm, the listed modes require a corresponding number of paths, which further increases the calculation burden although there may be countless possible combined modes in the entire multimodal network. From the user's perspective, they will only choose more competitive modes from a certain number of feasible alternatives. Tere is no need to list all modes for all OD pairs. Tese travel modes can be found indirectly from the path set. Te modal split can be built among the modes in the path set. Using this method, the number of variables can be constrained by the number of paths, which can signifcantly reduce the calculation burden. Following this idea, a generalized modal split/trafc assignment problem (GMS/TAP) is proposed in this paper.
Te second problem is to refect the balanced use of the entire network while cutting down the usage of the private car. Because the balanced use can be considered as no links and hubs becoming over congested, in this study, an exponential service-level function was formulated with the sum of the fow/capacity ratio of all the links and hubs. Meanwhile, in case the balanced use is achieved by increasing the car usage or other unexpected modal shift, a weighted modal shift value is formulated to guide the expected modal shift.
For the third problem, the IMNDP is a mixed-integer nonlinear bilevel NP-hard problem (MIBLP) [68]. Tere are usually two approaches to solve an MIBLP: exact search methods (exact optimization) and nonexact methods (heuristic), including genetic algorithms (GAs), evolutionary methods, ant colony optimization, tabu search, and simulated annealing [69]. However, in a real-size network, the IMNDP experiences a combinatorial explosion in the size of the solution space as the problem size increases. Te evaluation of potential candidate solutions requires the execution of a GMS/TAP, which requires a higher-order polynomial time [70]. To date, heuristic approaches have proved to be the only practical methods for solving realistic and large-sized problem instances [41].
Among the heuristic algorithms, the GA has been proven to be an efcient algorithm for solving an MIBLP because of its simplicity and minimal problem restrictions [71]. Many diferent types of GAs have been successfully adopted for MNDP [72][73][74][75]. Since IMNDP dealt with over three sets of decision variables, combining these diferent variables leads to challenges regarding infeasible solution management and neighbourhood evaluation. Tus, a hybrid genetic algorithm with adaptive diversity management (HGSADC) is proposed based on the work of Vidal et al. [76]. HGSADC is a hybrid metaheuristic that combines the exploration capabilities of GAs with efcient local search-based improvement procedures and diversity management mechanisms [77,78].
Te local search and diversity management in HGSADC is problem specifc and deals with vehicle routing problems with time windows. To adopt HGSADC for the IMNDP, a local search with a crossover, mutation, and repair operator; diversity management with a Hamming distance-based objective; and a dynamic penalty were introduced based on the previous research [73,79].
To summarize, this study proposed a new IMNDP concept to achieve a modal shift and balanced use of the entire network by the collaborative optimization of the hubs. Te main aim of this study was to provide new methodologies for the assessment and optimal design of a hub integration scheme. Te major contributions of this study are as follows: (1) Tis paper extended the MS/TAP to describe the mode choice for all modes in the path set, and a unifed optimization problem was formulated by combining the mode choice and route choice together. (2) A method for the optimal integration scheme was proposed, in which the modal shift and balanced use were chosen as the objectives, while the decision variables allow a collaborative optimization of hub locations, route layouts, and frequency setting for hubs with diferent modes. An HGSADC was proposed to address the optimal solution of the IMNDP.
Te abbreviations used in this paper are shown in Table 1.

Materials and Methods
In this section, the bilevel IMNDP is formulated for both the lower-level and upper-level problem. Te network representation is built in Section 3.1 at frst. Ten the formulation of lower-level GMS/TAP is divided into two parts: the modechoice structure and the user equilibrium (UE) principle. Te method of building the mode-choice structure is shown in Section 3.2. After combining the mode-choice structure with the user equilibrium principle, a unifed GMS/TAP is formulated in Section 3.3 along with its solution method in Section 3.4. Te upper-level model is introduced in Section 3.5, and the HGSADC and its operators are provided in Section 3.6.

Network Representation.
Five basic modes were considered in this study, i.e., bicycle, car hailing, car, bus, and metro. Te car, scheduled bus, and metro are the backbone of the urban transportation, complemented by bicycle and the online car hailing services such as Uber and Didi. Te mode set is denoted by M � bic, ch, car, bus, mtr { }. Denote the subnetworks as G m (N m , A sub m , H m ), where each subnetwork corresponds to a basic travel mode, m, where m ∈ M. N m , A sub m , and H m are the node, subnetwork link, and hub set, respectively. Tese subnetworks are combined with three link subsets. Te set of auxiliary links, A aux , connects the subnetwork nodes and their hubs within an acceptable walking distance. Te sets of embarking and alighting links (A emb and A ali , respectively) connect the OD nodes with the subnetwork nodes and represent the embarking and alighting activities. Te set of transfer links, A tra , connects diferent hubs and represents the transfer activities from one mode to another.
Te cost of the subnetwork link is the in-vehicle time, which is a link-fow-dependent function. Te alighting/ embarking link and auxiliary link comprise the walking time, waiting time, and time consumed in the hub, which is related to the length of the link, route frequency, and hub fow. Te impedance of the transfer link adds a pure penalty based on the above three times, which represents the user's natural resistance to transfer. Te impedance functions of the four links are continuously diferentiable and monotonically increasing functions.

Mode-Choice Structure.
To build a GMS/TAP, a modechoice structure is frst required to represent the travel modes. A modifed NL model is used in this paper, see Figure 1. Suppose the maximum transfer number is K. In level −1, all modes are classifed according to their number of transfers and divided into K + 1 columns (from 0 to K transfers). Ten, for each column k, a sub-NL model is built by modifying the NL structure for modes with k transfers.
Ten, from level 0 to level k, each basic mode of the k transfers combined mode is listed level by level. In the modechoice structure, nests in column k, level k are called mode nests (nests in blue boxes), which represent complete combined modes, while the nests in level 0 to k − 1 are called intermediate nests (nests in yellow boxes), which are the frst n modes of one complete combined mode.
To build a NL choice structure from the path set. A given path set of OD pair w is denoted as P w , where P w can be divided into K + 1 subsets P w k according to the number of transfers, k � 0 ∼ K. For the s th path p w k,s in P w k , s � 0 ∼ S w k . Te mode of p w k,s can be denoted with a mode chain mc w k,s (m 0 , m 1 , . . . , m n , . . . , m k ), where m n is the n th basic mode, n � 0 ∼ k. For convenience, a complete mode chain, mc w k,s (m 0 , m 1 , . . . , m n , . . . , m k ) is denoted as mc w k,s (k), and an intermediate mode chain with the frst n modes, mc w k,s (m 0 , m 1 , ..., m n ) is denoted as mc w k,s (n). For the NL mode-choice structure, denote the q th nest in column k and level n as e w k,q (m 0 , m 1 , . . . , m n ), which can be simplifed to e w k,q (n). In this study, if the mode sequence in mc w k,s (n) is the same as the mode sequence in e w k,q (n), then mc w k,s (n) � e w k,q (n). If the frst n modes in e w k,q (n + 1) are the same as e w k,q (n), then e w k,q (n + 1) ∈ e w k,q (n), meaning e w k,q (n + 1) is the lower-level nest of e w k,q (n). Following this corresponding relationship from mode chain to nest, for a given path set, the NL choice-structure can be built by matching the mode chains and the nests one by one.

Te Tree Conditions.
After building the mode-choice structure, the NL rules should be formulated among the fow of each nest. Meanwhile, to combine the mode choice and route choice together, this paper assumes that the path fow follows the user equilibrium principle. Tus, the following 3 conditions should be satisfed at the same time for a unifed formulation.
C1: Choice of Transfer Number. After the proportions are satisfed, no user will change the transfer number. Suppose the demand of OD pair w, g w is fxed. In the mode-choice structure, the fow of nests in column k, level −1 is denoted as g w k,q (−1). Ten, C1 is formulated as ∀w ∈ W, k � 0 ∼ K, and g w k,q (−1), and g w should satisfy the following constraints: where U w k,q (−1) is the equilibrium utility of the nest in column k, level −1, α −1 , and β −1 are the logit model parameters associated with level −1. C2: Choice of Mode. After the proportions are satisfed, no user will change the mode sequence. Denote the fow of each intermediate nest as g w k,q (n) and its lower-level nest fow as g w k,q (n + 1). C2 is represented as ∀w ∈ W, Journal of Advanced Transportation k � 0∼K and n � 0∼k − 1, and the following constraints should be obeyed: where U w k,q (n + 1) is the equilibrium utility of nest e w k,q (n + 1), U w k,q (n + 1) is the equilibrium utility set of all nests in level n + 1, and j w k,q (n + 1) is a binary variable for determining whether e w k,q (n + 1) ∈ e w k,q (n).
Equation (4) states that only when path cost c w k,s equals its corresponding mode nest's equilibrium utility U w k,q (k), then can path p w k,s carry fow; otherwise, the path fow is 0. It should be noted that the condition mc w k,s (k) � e w k,q (k) is to fnd the corresponding mode nest e w k,q (k) of the complete mode chain mc w k,s (k), and only paths in the same mode nest should obey the UE principle.
However, the variable of C1 and C2 is the nest fow, while for C3, it is the path fow. To keep the variables consistent in the three conditions so that a unifed formulation can be built to combine the mode choice and route choice together, a simple transformation is built.
Denote g w k,s (−1) as the fow of path p w k,s 's corresponding nest in level −1 and g w k,s (n) as the fow of mc w k,s (n)'s corresponding nest in level n, which means g w k,s (−1) � g w k,q (−1) and g w k,s (n) � g w k,q (n). Ten, C1 can be reformulated as While C2 can be reformulated as ∀w ∈ W, k � 0∼K, s � 1 ∼ S w k , n � 0∼k -1, and mc w k,s (n) � e w k,q (n): Equations (5) and (6) only replaces g w k,q (−1) and g w k,q (n) in equations (1) and (2) with g w k,s (−1) and g w k,s (n), and the other variables are the same.

Unifed Optimization Problem.
Denote the equilibrium fow vector of all OD pairs as F * . Te following variational inequality (VI) formulation was built to unify the above three conditions. Te proof is presented in Appendix A.

Proposition 1. Te optimal solution F * of the proposed VI problem gives the NL mode choice solution and UE route choice solution together only if there exists the set of values λ
where Te solution space of the VI problem is Equation (9) ensures that the total path fow of OD w equals demand g w , and equation (10) indicates that the fow is nonnegative.
Te VI problem can be reformulated as a single-level NLP according to Aghassi et al. [80]. Since the solution space (equations (9) and (10)) of the VI problem can be reformulated as a nonempty polyhedron, . Te VI problem in equation (7) can be transformed into the following NLP.

Proposition 2. Suppose that D is a nonempty polyhedron,
Ten, F * solves the VI problem in equation (7) if the following optimization problem has an optimal value of 0 and ∃ξ * ∈ R W , such that (F * , ξ * ) is an optimal solution: min Journal of Advanced Transportation

Hybrid MSA Algorithm.
Te NLP can be solved using the currently available software such as CPLEX, Gurobi, and MATLAB fmincon solver. Meanwhile, the standard formulation provides an excellent coverage criterion, where the minimum objective function should be zero. However, the number of variables in the two vectors, F and ξ, which are related to the number of paths and number of OD pairs, respectively, can still be very large in real-size network. A highly efcient algorithm is necessary to apply the GMS/ TAP to the IMNDP. Tus, a method of successive averages (MSA)-based method is proposed, which is commonly used to solve trafc assignment problems [81].

Formulation of the Upper-Level Problem.
Te notations used in the upper-level problem are shown in Table 2. IMNDP aims to realize a modal shift and obtain a balanced use among the entire multimodal network under budget limits. Tis results in three objectives: construction cost, modal shift, and balanced use. Because the balanced use can be considered as no links and hubs becoming over congested, an exponential service-level function was formulated with the sum of the fow/capacity ratio of all the links and hubs. Meanwhile, in case the balanced use is achieved by increasing the car usage or other unexpected modal shift, a weighted modal shift value is formulated. As a widely used objective, the total time of all users was also adopted to protect the users' interests. Ten, the objective function of the problem is Since it is beyond the work of this paper to cover all hubs of diferent modes, for simplifcation, only the fve modes in Section 3.1 are considered. Te optimization of hubs is classifed into two types. One is the location-based optimization, which are parking lot and bicycle station. Te optimization only needs to determine their locations. Te other is the route-based optimization, which is the bus stop. Te change of their hub locations will also infuence the route layouts and frequencies. Te integration scheme is to optimize these variables simultaneously to build a wellconnected multimodal network. Te objective functions are calculated with equations (13)-(16) using these results: Equation (13) is the total travel time. Equation (14) is the total construction cost, which consists of the demolition and construction costs of the hubs, the bus purchase cost, and the bus operational cost. Equation (15) is the modal shift value. To analyse the modal shift of OD pairs with diferent lengths, mf l,m is denoted as the total fow of mode m among all OD pairs with distance l. Te predefned weight parameter ms l,m is used to value the 1 unit of fow change compared to the original network if ms l,m is set as a minus value, which means the modal shift of this mode in this distance is encouraged to increase, and vice versa. Equation (16) is the service-level value, which is an exponential service-level function with the sum of the fow/capacity ratio of all the links and hubs. To ensure that highly loaded hubs or links will have higher service level value penalty, a parameter θ is added in equation (16).
For location-based optimization, the constraints include those as follows: Equation (17) ensures that for each candidate node, only one hub with mode m can be placed. Equation (18) is the hub number constraint for hubs with mode m. 8 Journal of Advanced Transportation For route-based optimization, where the decision variables are bus layouts and frequencies, referring to [8,73], the constraints include those as follows: Among the constraints, equations (19) and (20) ensure that all bus routes start and end at a bus terminal. Among OD pairs with direct distance between l − 1 and l km, the total fow of mode m in the optimization scheme or in the original network sl h m , sl a m Te service level value of the hub or link with mode m X ojr 1 if route r starts at node j, and 0 otherwise X ior 1 if route r ends at node i, and 0 otherwise X ijr 1 if route r (r � 1toR max ) passes through node j immediately after node i, and 0 otherwise  (21) ensures that each bus interchange has one preceding node and one subsequent node. Equations (22)-(24) ensure that each bus node can be visited by one particular bus route only once at most (bus node repeat constraint). Equations (25)- (28) are the bus stop distance, bus route length, bus frequency, and bus feet size constraints, respectively. Te bus frequency of one route is calculated using the assigned number of buses in equation (29).

Solution Method of the Upper-Level Problem.
Te upperlevel problem is a MIBLP with 4 sets of variables. To solve the problem, the HGSADC framework is introduced as follows.
Te local search and diversity management in HGSADC is problem specifc. Tus, the local search with a crossover, mutation, and repair operator; diversity management with a Hamming distance-based objective; and a dynamic penalty were introduced one by one in the following sections.

Initialization Operator.
For route-based optimization, the initialization operator determines the initialization locations of the bus stops and their sequences in each bus route. Te initialization scheme should satisfy the stop distance constraint, route length constraint, and bus node repeat constraint.
Te bus route frequency is initialized by randomly assigning buses to each route according to equation (29). If the frequency violates frequency constraint and the number of buses violates feet size constraint, the process will repeat until the constraints are satisfed. It should be noted that lines 4 ∼ 12 is called the bus sequence generation function in the rest of paper. Tis function can be used to generate a bus sequence with defned starting and ending nodes.
For location-based optimization, the parking lot and bicycle station initialization should determine the number and position of each hub. Te initialization operator randomly adds or deletes several hubs and ensures that the total number of hubs satisfes the hub number constraint.

Crossover
Operator. Te crossover of the bus route consists of two levels: combinations of diferent routes and diferent sequences of bus stops in one route. To increase the exploration capability and ensure the satisfaction of the constraints for the generated routes after crossover, referring to Ngamchai and Lovell [79], two crossover methods for route-based optimization, i.e., line crossover and stop crossover, are used. Line crossover exchanges two diferent routes in two solutions, while the stop crossover operator is triggered when two routes in one solution have segments that can be spliced.
For location-based crossover, two hubs in two solutions are randomly selected; if they are diferent and not repeated in the other solution, then exchange.

Mutation Operator.
Te mutation process generates diferent genetic materials to facilitate exploration in the search space. According to Szeto and Wu [73], four types of mutation operators for route-based optimization are proposed: insert, delete, swap, and transfer.
Te insert mutation is used to insert an interchange node in a randomly selected route. Te deletion mutation is similar. Te swap mutation is to exchange the nodes of two routes if the stop distance constraint can be satisfed after the swap, while transfer mutation is to transfer one node to the other route, and similarly, the stop distance constraint should be obeyed after the transfer.
Te location-based mutation only contains the insert and delete mutations. Te insert mutation inserts a random hub that is not in the solution, whereas the delete mutation deletes one random hub in the solution.
3.6.4. Repair Operator. After the crossover and mutation processes, the generated children may violate the constraints, and a repair operator is needed to fx these infeasible solutions. Te constraints for route-based optimization are route length, stop distance, node repeat, and frequency, while the constraint for location-based optimization is the number of hubs. With four constraints, there is a high probability that the bus route will violate other constraints after one constraint is repaired. A new four-stage repair strategy was proposed in this study. Te solution set will undergo four repair operators in the sequence of the route length operator, stop distance operator, node repeat operator, and route frequency operator.
(1) Route Length Repair Operator. See Figure 2, when route r 1 is too short, the algorithm selects two nodes, n 1a and n 1b , with a sufcient distance and inserts an intermediate node n 1c between them. Ten, two stop sequences from n 1a to n 1c and from n 1c to n 1b are generated with the bus stop sequence generation function in Algorithm1, and the original stop sequence between n 1a and n 1b is replaced. Te situation when the route is too long is similar, and two nodes with a sufcient distance are chosen and replaced with the generated stop sequence using the bus stop sequence generation function. (2) Stop Distance Repair Operator. A stop sequence that violates the stop distance constraint may contain only two nodes or consist of several consecutive nodes in which all the adjacent nodes violate the stop distance constraint. When the stop distance is too long, and the stop sequence has only two nodes, see Figure 3, the algorithm fnds an intermediate node to insert between the two nodes. When the stop sequence has over two nodes, a new stop sequence will be generated to replace it. Te situation when stop distance is too short is similar. (3) Route Repeat Repair Operator. For repeated nodes, one node is randomly selected to be unchanging. Ten, for each repeated node, n 1a , a node is selected from CN that satisfes the stop distance constraint with both n 1a − 1 and n 1a + 1 to replace n 1a .
(4) Route Frequency Repair Operator. Te frequency repair operator adds or deletes one bus until the bus frequency constraint is satisfed.

Diversifcation Mechanism and Infeasible Solution
Management. To ensure the diversity of the solutions, the Hamming distance is introduced to represent the diferences between one solution and the best solution, which is also called the diversity contribution function ΔS. It is defned as the average distance from S to its closest neighbours in the subpopulation. ΔS consists of three parts, which are the bus route, parking lot, and bicycle station Hamming distance. According to Szeto and Wu [73], the bus Hamming distance is the number of diferent consecutive node pairs for each pair of corresponding routes compared to the best solution. For the location-based one, it is the number of diferent hubs compared to the best solution, and ΔS is the sum of the above three distances. A dynamic penalty is then introduced to handle the infeasible solutions. Suppose the number of feasible and infeasible solutions are N inf and N fea , respectively. A basic penalty value, ϖ, is adjusted for every μ iterations, and when the number of N inf increases by 1, ϖ decreases by λ%, and vice versa. (see Algorithm 2) (1) for route r from 1 to R max (2) Whether route r needs to construct a new terminal with probability P s ter and P e ter (3) Determine the starting and ending node: if new terminal should be built, select a node from CN as the starting node n s r or ending node n e r of route r; if no, then select a node from TN as the starting node or ending node (4) while the t th node of roue r, n r (t) is not the n e r , and the route satisfes the route length constraint, bus node repeat constraint, and time < T max do (5) Find the node set N stop r (t) with distance to n r (t) satisfying the stop distance constraint (6) if N stop r (t) is empty or all nodes have been visited in N stop Randomly select a node n stop r (t) from N stop r (t) that has not been visited (10) if |n e r − n temp   Ten, the survivors of each iteration are selected based on the rank of objective function fit(S) and the rank of Hamming distance h dr(S), which is called the biased ftness in this paper. For infeasible solutions, the objective function should multiply dynamic penalty values, ϖ. Referring to Vidal et al. [76], the biased ftness, BF(S), is formulated: where N elite is a parameter ensuring elitism properties during survivor selection. After building the whole bilevel model, the extracts of the algorithms can be found in Appendix B.

Results and Discussion
Tis section discusses three experiments that were performed to validate the proposed models and algorithms. In experiment 1, the exact results of the standard NLP of the GMS/TAP were verifed, along with the nonexact results of the hybrid MSA algorithm. In experiment 2, the results of the upper-level problems were validated by testing the diversifcation mechanism and biased ftness function. Te application of the IMNDP in a real-size network was then tested in experiment 3.

Experiment 1. Te frst experiment was conducted to
show the solution properties of the GMS/TAP. Te ninenode network used in Wu et al. [82] was adopted, as shown in Figure 4. Te link impedance functions of the road links and bus links were BPR (Bureau of public road) functions. Te link free-fow travel time and capacity were noted in the two-tuple next to each link.
To adopt the original network into GMS/TAP, a modifed network with the same properties was proposed. As shown in Figure 5, numbers 1-8 are the road nodes, 9-17 are the bus nodes, and 18-21 are the OD nodes. Te transfer links and embarking/alighting links were added at a cost of zero. Te road and bus subnetwork links were maintained with the same parameters. Meanwhile, although the car + bus combined mode is difcult to observe in real life, they were allowable in this experiment, and car users could directly transfer from a car node to a bus stop, and vice versa.
For simplifcation, the car and bus modes are denoted as c and b. Although it was unnecessary for a path set to cover all possible modes, in this experiment, six modes were selected to show the solution properties of the GMS/TAP. Te corresponding modes included those as follows: (1) 0 transfers, where the travel mode is either c or b; (2) 1 transfer, where the travel modes are either b + c or c + b, and (3) 2 transfers, where the travel modes comprise b + c + b and c + b + c. Tree paths for each mode with the lowest path costs were selected. If there were fewer than three paths, only the available paths were selected.
Te model was calculated when maximum transfer number K � 0, 1, and 2 using the MATLAB fmincon solver. Meanwhile, when K � 2, the results were also calculated using the hybrid MSA algorithm. Te logit parameters at each level were α −1 � 0.5, β 1 � 0.85, α 2 � 0.5, and β 2 � 0.8. Te path set for is presented in Table 3. A bold number in the path column represents the transfer hub. Te optimal fow, F * , and optimal objective, ξ * , makes the NLP have the optimal objective, and equilibrium path cost when K � 2 is based on the fmincon solver, and the hybrid MSA are listed in Tables 4 and 5.
In Table 4, the optimal objectives based on the fmincon solver is 8.6484e − 6, which is close to the coverage criterion zero of Proposition 2. Meanwhile, for all the paths with the same mode, such as paths 18,19, and 20 with mode c, the paths carrying fow all have the lowest path cost, so the UE principle is satisfed.
Te other two conditions in section 3.3.1 are also verifed. Taking OD 18-20 as an example, the minimum utility values of 0, 1, and 2 transfers are 55.0, 59.6, and 58.0, respectively, and the equilibrium fows are 28.7, 0.1, and 1.2, which obeys condition 1. Among the two transfer trips, the minimum utility values are 59.6 and 58 for modes b + c + b and c + b + c, respectively, and the total fow of 1.2 is divided into 0.2 and 1.0, which also obeys condition 2. Table 5 shows that the objective function is 3.5948, and the UE principle is also obeyed. Compared with       Figure 6. Te three subnetworks from left to right are the road, metro, and bus subnetworks. Te link and hub free-fow travel times and capacities are denoted by blue and green texts. Te following modifcations were made based on the network of experiment 1: (1) A new parking lot was added in node 26 and connected to road node 5 with bidirectional auxiliary links 37 and 51. Users must transfer from a parking lot to a metro station to execute a P&R combined trip. (2) Two parallel metro lines with a frequency of 0.083 h/ veh were added to provide more mode choices. Only one bus route 21-25 was kept at 0.068 h/veh to leave space for optimization.
(3) In experiment 1, the costs of transfer links, embarking/alighting links, and auxiliary links were default to be 0. In this experiment, to verify the results after the optimization of the hubs, the costs of transfer links, embarking/alighting links, and auxiliary links were added with the sum of the walking time of the links, walking time inside the hubs (BPR function with green text as the parameters), waiting time (the route frequency), and pure penalty (a fxed time of 0.083 h).
(4) A Euclidean space coordinate system was built, where the x and y axes represented the position of the node (see Figure 7). Because the four OD pairs were connected to several direct car, bus, or metro services to avoid the situation where no users take combined trips, one unit of length in the network was set at 10 km, and the car, metro, and bus speeds were set at 30, 25, and 20 km/h, respectively.
Experiment 2 aimed to optimize the parking lots and bus routes to integrate the network. Te maximum number of bus routes was set to two, while the number of parking lots was set to 1-6. Te original number of buses was fve. Te weights of the four objectives such as the total time, construction cost, modal shift value, and service level value, were set at 1, 1, 2, and 3, respectively.

Results of Experiment 2.
Four sets of GAs tested the model performance: GA1 was the normal GA, GA2 was a GA with a diversifcation mechanism alone, GA3 was a GA with a biased ftness function alone, and GA4 was the HGSADC. Each GA was ended when 500 generations were produced and repeated for 20 runs. Te average results of the original network and 4 GAs are listed in Table 6.
In all 20 runs, GA2, GA3, and HGSADC all converged to the same results with a weighted objective value of 17.3641. However, GA1 was trapped in a local optimal objective with an average objective value of 21.1058, showing that both the diversifcation mechanism and biased ftness function contributed to improving the exploration capacity. Meanwhile, HGSADC had far fewer convergence iterations than the other three GAs. Because the diversifcation mechanism involved a reinitialization process to obtain more genetic materials, it also made the average CPU time longer than the other three GAs. Considering it is the total time of 500 iterations, the fast convergence speed can be useful in large-scale problems.

Analyzation of Schemes.
Te performances of the original, GA1, and HGSADC schemes were then compared. See Table 7.
(1) Original Network. Te fow distribution of the original network is shown in Figure 8. From left to right, the multimodal network comprises a green road subnetwork, brown metro subnetwork, and a purple bus subnetwork. Te link fow is noted beside the link in black text, while the hub fow is shown in blue text. Meanwhile, transfer fow is noted beside the transfer link. Te pink text shows the mtr + car and mtr + bus transfer fow, while the orange text shows the car + mtr and bus + mtr transfer fow. Te node fow is the sum of all the links entering the node. Te two parallel OD pairs, 9-11 and 10-12, can be served by the car and two parallel metro lines directly. Te speed of the car is faster than that of the metro, while the capacity of the road link is 10, which is not sufcient to meet the 30 demands of the parallel OD pairs. Tis portion of the users will choose the single-car mode frst and then the metro. Meanwhile, because a parking lot is built at node 26, and road link 1-5-7-3 is congested (the fow/capacity ratio is higher than 1), while metro line 13-14-15-16 is underutilized (lower than 0.75), some users will transfer from parking lot 26 to metro station 14 to bypass the congested road link, 5-7-3.

Journal of Advanced Transportation
For the two diagonal OD pairs, OD pair 9-12 can be served by the bus route directly. Among the 30 demands, 10.2 users choose bus link 21-22, and the other 8.9 fow of this route is the mtr + bus fow. Te bus route relieves the congestion of the road links in this direction, and road links 1-6, 5-6, 7-8, and 7-4 all have a good service level. However, the other diagonal OD pair, 10-11, can only be served by the road subnetwork directly, which makes the service level of related links 2-5 and 8-3 very low.
(2) GA1 Scheme. Te GA1 scheme adds a new parking lot in node 31, as well as bus route 26-27-2829 with frequency 0.12 h/ veh. Te original parking lot 30 and bus route 21-22-23-24-25 are kept, but the frequency of the original bus route increases to 0.112 h/veh. Te fow distribution is illustrated in Figure 9.
Te newly added bus route, 26-29, only attracts 3.7 users because it is parallel to metro route 13-16. However, the congestion of road links 1-5-7-3 is still slightly reduced. Meanwhile, the new parking lot at 31 attracts 14.6 users to use the P&R mode, which reduces the congestion of road link 6-8-4.
Te added bus route provides a direct bus service for OD pair 10-11, which attracts eight single bus mode users to bus stop 26. Tis bus route creates two new mode choices for this OD pair, mtr + bus and bus + mtr, and except for 8 direct users, 8.1 users transfer from metro station 18 to bus stop 27 and 2.2 users transfer from bus stop 29 to metro station 13, which signifcantly relieves the congestion of road links 2-5 and 8-3. With more users choosing PT, the fows of other single-car mode paths, such as 2-5-7-3 and 2-6-8-3, also decrease by a certain amount, which leads to an overall reduction in congestion of road subnetwork. Because the weight of the service level value is the highest, the lowest service level value of 405.4 makes the weighted objective of the HGSADC lower than that of the GA1 scheme even though the total time is similar, and the construction cost is higher.   Since Nanjing is divided into two parts by the Yangtze River, Linjiang Station is an important node which connects the south and north banks of the city, and a typical P&R parking lot is built there. Te network comprises 5 subnetworks, i.e., bicycle, car, car hailing, metro, and bus subnetwork. In Figure 12, from (a) to (e), the blue, green, and red subnetworks are bicycle, car, and car hailing subnetworks, and they share the same physical network with 125 nodes and 428 links; meanwhile, car and car hailing subnetwork share the same fow. Te blue bicycle subnetwork contains 55 bicycle stations while the green car subnetwork has 45 parking lots. Te brown metro subnetwork has 2 metro lines and 18 metro stations. Te purple bus subnetwork consists of 11 bus lines and 95 bus stops. Te district is divided into 22 trafc zones with 462 OD pairs, see Figure 12(f ).

Results of Experiment 3.
Te decision variables were bicycle stations, parking lots, bus layouts, and their frequencies. Te iteration number was set as 3000, and the maximum transfer number K was 2. Te number of bus routes was kept at 11. Te number of bicycle stations, parking lots, and bus stops was set from 50 to 200, 45 to 100, and 80 to 200, respectively. Te bus stop distance constraint was from 0.3 km to 2 km and 0.0833 h/veh to 0.5 h/veh for bus route frequency constraint. Te bus route length constraint was from 5 km to 25 km. 40 paths for each OD pair were chosen to ensure the stability of hybrid MSA in obtaining the equilibrium fow. Te pure penalty is 5 mins.
(1) Travel Time and Construction Cost. Te algorithm converges in 2335 iterations with 172417.54 s. GMS/TAP function accounts for 77% of the total time, in which path generation, path utility, and hybrid MSA takes 39.5%, 4.8%, and 32.7%. Te HGSADC and other functions take 23% of the total time. Te best solution contains 51 parking lots, 74 bicycle stations, and 11 bus routes with 117 bus stops. See Figure 13. With a construction cost of 11566314.32 RMB, the travel time decreases from 749042.8 h to 680345.5 h, which improves 9.2%. Since this study concentrates on the modal shift and balanced use of the whole network, comparing with Szeto et al. [69], in which a 22.7% reduction of travel time in transit subnetwork is obtained, the improvement of travel time is not obvious in this experiment. 17       (2) Service Level Value. Te service level of the original network is presented in Figure 14, and in the original network, the road and metro subnetworks are overloaded, especially in the CBD area with longitude from 118.7 to 118.72, latitude from 31.98 to 31.99. Meanwhile, the fow between north bank and south bank also concentrates in the car and metro mode, and the service level in the corresponding links and hubs are all lower than level E. With a congested metro and car subnetwork, the bicycle and bus subnetworks in Figures 14(a) and 14(e) are underutilized; to fully use the potential capacity of these subnetworks, a modal shift is necessary.
In the optimal scheme (see Figure 15), the service level value decreases from 40253.82 to 29903.51, which improves 25.71%. Although the main trafc corridor from west to east is still congested in road subnetwork, the number of over congested links and hubs (service level higher than 1) decreases with 21 and 3 for road and metro subnetworks, respectively. Tese fows transfer to bicycle and bus subnetworks, and an increase in usage with 36 and 44 links and hubs is obtained. In the x axis, the mode bic, car, ch, mtr, and bus are denoted as b, c, h, m, and B, respectively. Te combined trip, such as car + mtr is denoted as cm. While on the y axis, the mode fow is classifed per kilometers according to the direct OD distance. Te z axis is the amount of mode fow.
In the original network, 77.33% of the fow is singlemode within which 37.17% comes from the private car mode. Since there are only two metro lines while a complete bus network is built, the metro subnetwork is mostly used in combined trips. Te bicycle is normally used alone in short distance trips and in long-distance trips to complement bus and metro subnetworks, and the car hailing subnetwork plays a similar role as the private car but with a lower fow.
Considering the original road subnetwork is over congested, the single-mode c is encouraged to be reduced with positive weight value. Meanwhile, one transfer modes in short distance and 2 transfers modes in short and middle distance needs to be cut down as well. Te modal shift weight matrix is divided into seven columns and three rows, see Table 8.
Te fow of each mode in the optimal scheme is shown in Figure 17. Te 0 transfer modes and two transfers modes 32      users decrease with 6885 and 12797, respectively, and these users transfer to one transfer mode. Among single-mode trips, a decrease of 27768 car users and an increase of 4853, 14860, and 20853 bic, mtr, and bus users are observed. Tese four single-modes produce a weighted modal shift value −163681. In one transfer trips, the modal shift does not change as wishes, the increase in 19682 one transfer modes users concentrate in short distance OD pairs, especially the bm and mb modes. But in other perspective, this means the new bicycle station scheme improves the connection between bicycle and metro subnetworks. For two transfers modes, the low fow does not have an obvious infuence on the modal shift value.

Sensitive Tests of the Parameters
(1) Efect of the Transfer Penalty. Tis section tested the efect of transfer penalty on the fow distribution of the optimized network. Since there are several types of transfer costs which are related to the connected modes, it is difcult to test these parameters one by one. Considering every transfer involves in a pure penalty with 5 mins, pure penalty represents the natural resistance of the users to the transfer. Tus, the experiment tested the network with a pure penalty from 0 to 15 min, with an interval of 0.25 min. Te total mode miles are chosen as the index, which is used to describe the overall usage of the corresponding subnetwork. See Figure 18.
Journal of Advanced Transportation In Figure 18(a), when the transfer penalty increases, the miles on car and bicycle subnetwork increases. For the metro and bus subnetworks, which rely more on transfer, the increase of transfer penalty leads to a decrease of their travel miles. In Figure 18(b), we can see that the number of transfers is highly related to the transfer time. Intuitively, the decrease of the number of transfers is good for the transportation system. However, in this test, this is achieved by increasing the transfer time rather than improving the connection of the network. Te long-time of waiting pushes the users to the congested links of the single modes, while they cannot use the underutilized links of other subnetworks through combined trip to bypass the congested links. Te results suggest that the number of transfers alone may not be a good index to value the integration of the multimodal network.
(2) Efect of the Congestion. We further investigated the congestion efect under diferent demand levels. See Figure 19. Te OD demand varies from 5% to 300% of the original demand with an increase of 5%. As the OD demand increases, the users choose road and metro subnetworks at frst. When the demand reaches 45% of the original demand, the car subnetwork becomes relatively more congested than other three subnetworks, i.e., bus, bicycle, and car hailing, and their share of travel miles increases while car and metro decreases. As shown in Figure 15, under the original demand, the car and metro subnetworks are already over congested, the decrease speed for these two modes is increasing, and more users choose other three modes.
Tis result suggests that there may be one or two modes that the users will choose at frst because their competitions are   Step 1 Step 2 Step 3 Step 4 Step 5 Step 6 much better than other modes. Only when these preferred subnetworks are congested, then users will choose the remaining modes. Te natural diferences of the travel modes are hard to change only via an improvement of the infrastructures. A systematic demand control and operational management are also required to balance the fow distribution.

Conclusions
Tis paper proposed a bilevel IMNDP to optimize subnetwork hubs, guide a modal shift, and achieve a balanced fow distribution among all the subnetworks. Te upper-level problem had the goal of optimizing various network hubs, and the decision variables were classifed as the route-based and location-based optimization. Te route-based optimization included the route layout and frequency setting, whereas the location-based consisted of the hub locations. In the lower-level problem, a NL mode-choice structure was devised to describe the mode choice of the contained modes in the path set. A unifed VI formulation was developed and transformed into an NLP that could be solved using commercial software. To facilitate the application in a real-size network, a hybrid MSA was further developed to solve the lower-level problem. To solve the proposed bilevel model, an HGSADC along with the initialization, crossover, mutation, repair operators, biased ftness, and diversifcation mechanisms were introduced in this paper. Te solution properties and model performance were verifed in three experiments, and a signifcant modal shift and service level increase could be observed. Tis study opens several research directions. (1) In this study, the modal shift was valued with a predefned matrix, which was based on the current use of the network and experience, while the balanced use of the network was assessed by the service levels of all the links and hubs. However, how to assess a real balanced use and what is a reasonable modal shift pattern still requires much work. (2) Te optimization of different hubs in this paper only considered limited number of constraints. However, the metro, car hailing, and many other modes were not shown in this paper, and the optimization of these modes and their combinations also requires deeper research. (3) As discussed in the literature review, the integration of multimodal network can be achieved by many other factors, such as policy, environment, market needs, fare, and coordinated schedules which can be used as the decision variables as well. (4) Te impedance function requires a systematic data collection and calibration methods regarding to the diferent links, hubs, single, and combined modes. (5) Besides the mode choice, multiple classes of users, link interactions, dynamic assignments, demand uncertainties, exact-method algorithms and many other factors could also apply to the GMS/TAP.
Te main contributions of this paper are as follows: (1) the extension of the traditional MS/TAP model framework for a wider range of modes. (2) Te methodologies for the assessment and optimal design of a hub integration scheme. However, these contributions are more theoretical than practical. Without a calibrated impedance function, the model is not suitable to be adopted into the real-size network with more travel modes. Tus, the priority of our future research is to calibrate a system of impedance functions. Considering there are mature data for the car, metro, and bus subnetworks separately, the research could start from the extraction method of the transfer fow between the three modes using the cellular signaling data, bus card data, and GPS data.
Step 4: GMS/TAP Process. Each solution in the feasible and infeasible subset should go through a GMS/TAP process. Obtain the NL mode-choice structure of each OD pair, then call the Algorithm 3 to output the objective functions of each solution.
Step 5: Biased Fitness. Calculate the Hamming distance of feasible and infeasible subset, and obtain the biased ftness by adding the hamming distance and objective function. Select the survivors and set the iteration number t � t + 1; Step 6: Diversifcation Mechanism. Determine whether the algorithm should end and whether a diversifcation should be conducted.

Data Availability
Te data can be achieved by contacting the frst author.

Conflicts of Interest
Te authors declare that they have no conficts of interest.