Multicriteria Intermodal Freight Network Optimal Problem with Heterogeneous Preferences under Belt and Road Initiative

: In this study, we demonstrate the importance of incorporating shippers’ preference heterogeneity into the optimization of the China Railway express network. In particular, a bilevel programming model is established to minimize the total construction cost for the government in the upper level and maximize the shippers’ satisfaction in the lower level. The proposed model considers price, time, reliability, frequency, safety, ﬂexibility, traceability, and emission. Two designs are obtained by applying the model to two scenarios, in which one is of the aggregate shipper group and the other is of the three distinct clusters. Results show that explicitly including heterogeneity in network optimization pays o ﬀ in terms of the dramatic increase in shippers’ satisfaction and the share of the sustainable railway without generating extra cost for the system. The results of this study could lead to insightful implication for proper network planning for the China Railway express and some useful suggestions on the subsidies of the government. This study aims to explore the decisions regarding achieving an e ﬃ cient network of CRexpress with mass customization service. When designing the CRexpress service network, the LSI determines an optimal conﬁguration in terms of terminal locations and allocates the shipper demand to di ﬀ erent LSPs.


Introduction
The Belt and Road Initiative (BRI), announced in 2013, together with the "Vision and Actions on Jointly Building the Silk Road Economic Belt and 21st-Century Maritime Silk Road" has provided many opportunities for the trade between China and Europe along with huge benefits for China Railway Express (CRexpress). Since then, China has faced many severe environmental issues because freight transportation is the main contributor of emission.
Over the past 8 years, the number of CRexpress lines experienced a rapid growth from 17 lines in 2011 to 6363 lines in 2018. To date, 56 cities operate CRexpress lines, reaching 49 cities in 15 countries in Europe. In 2015, the Leading Group Office on the Construction of the Belt and Road issued the Development Plan of China-Europe Freight Train Construction (2016-2020), which confirms the significance of the development of the CRexpress given that it serves as an important carrier and gripper to promote the construction of the "Belt and Road". Apart from the thrilling progress, the development of CRexpress has diverse challenges. One key bottleneck that has to be eased is how to establish an environmental sustainability-oriented coordination mechanism and optimize the overall layout of the CRexpress system logistics service network.
The BRI service network has three decision makers: the shippers, the railway companies, and the central government. The shippers aim to get their cargo shipped with abundant benefits. The central government, the logistics service integrator (LSI), will distribute the orders from the clients to the logistics service providers (LSPs), which are the local railway companies. The LSPs will fulfill the demand by providing the corresponding customized services. At present, LSPs operating the lines only aim to be at an advantage in the competition with other LSPs; therefore, services with quite low price are supplied to attract shippers. At the beginning, this approach helped expand the volume, but it also brings bad results as the lines increased. First, many lines covering the same area leads to a huge waste of resources. Second, the cargo supply is insufficient due to the vicious competition. Finally, as a result of the above two scenarios, the LSI has to give a large amount of subsidies to LSPs. Under such a circumstance, the LSI of the CRexpress system guides the development of the LSPs from the overall level and improves the operational quality.
To solve the problem faced by the CRexpress system, a well-researched concept is that of service network design problem (SNDP). The SNDP is a key typical tactical problem in multimodal transportation, which mainly focuses on how to provide proper logistics services. Luathep et al. [1] proposed a global optimization algorithm to solve the transportation network design problem, which is expressed as a mathematical programming with an equilibrium constraint. Alumur et al. [2] addressed a hub location problem with different traveling modes and service time promises. Ambrosino and Sciomachen [3] proposed a mixed integer linear programming (MILP) model in which the containerized demand flow can be split into different services to evaluate suitable sites for dry-ports in northwestern Italian regions. Munim and Haralambides [4] analyzed the cross-border competition and cooperation between ports with a mixed integer linear programming. Their work shows that port users would greatly benefit from the cooperation, although some port authorities would suffer a revenue loss, which could be somehow compensated. Yang et al. [5] dealt with the reconstruction problem of the shipping service network between Asia and Europe. In particular, the authors tested the proposed bilevel model under different scenarios to obtain the new optimal networks. Some researchers investigated the CRexpress problems. Jiang et al. [6] were one of the first researchers that stressed the benefit of setting up consolidation centers for CRexpress. They built a mixed integer programming model, which selected Urimqi to be the only hub for 34 cities all over China.
The above contributions to SNDP are made to achieve a system-optimal solution. However, modeling for freight transport demand has significant wields on SNDP with the expanding commercial competitiveness over the past decade [7]. During a demand analysis, many sophisticated disaggregated models based on individual data have been used. Travel time cost was considered the sole objective in early transportation planning, however, people realized that a transportation network is a large system confronted with varied needs from different aspects [8]. Vinod and Baumol [9] first developed a model with shippers' choice to describe the mode choice with attributes of cost, time, reliability, and safety. Since then, different attributes have been investigated to explain the problem related to freight transport choice. Researchers realized that simply allowing for factors such as cost and time is not enough [10,11]. Cullinane and Toy [12] found out that "cost", "time", and "reliability" are the three commonly used attributes. Shinghal and Fowkes [13] conducted a stated preference survey on the Delhi to Bombay corridor. Their results show that frequency is of great importance to the mode choice. Yu et al. [14] claimed that flexibility could have distinct effects on logistics services under different environmental conditions. Tuzkaya and Önüt [15] evaluated the modes between Turkey and Germany in which traceability is also selected. Transport has great responsibility for environmental impacts, thus, emission has gained great attention-many articles are concerned about emission [16]. Diverse literature with multiple criteria could be found in Table 1. On the basis of the analysis above, we selected eight attributes to reflect the shippers' demand preference: cost, time, reliability, frequency, flexibility, traceability, emission, and safety.
Methods, such as multicriteria decision-making (MCDM), can be used to unify more than two indices, and approaches of different types of MCDM are employed. Liu et al. [17] used the technique for order preference by similarity to ideal solution (TOPSIS) model combined with entropy weight to evaluate the nodes' importance and features of Shanxi water network and Beijing subway network. Long and Grasman [18] provided a multicriteria decision framework to evaluate the location of inland freight hubs, in which the relevant criteria needed were selected from subject matter experts. Zhao et al. [19] proposed a relative-entropy rank evaluation method for multiple-attribute decision-making. Hui et al. [20] also used multicriteria decision-making in synthesizing the multiple indices to evaluate a complex network. Wang and Yeo [21] built an evaluation structure for intermodal routing from Korea to Central Asia by adopting the Fuzzy Delphi and Fuzzy (Elimination and Choice Translating Reality I) methods. Their cases include five principle factors, among which total cost turns out to be the most important factor that affects the companies when they are selecting a route. We also summarize some applications of integrated attributes and shipper preferences under BRI. Jiang et al. [22] examined the choice probability of two types of goods between CRexpress and shipping among five typical CRexpress routes with a binary logit model. They conclude that government subsidies greatly contribute to decreasing the cost of CRexpress operators, and logistics service producers are more willing to choose CRexpress than others. As an extension of this work, Zhao et al. [23] evaluated cities from the perspective of government policy and CRexpress operation experience with TOPSIS and used mixed integer programming to select the key sites for consolidation. Studies that considered multicriteria decision-making are summarized in Table 1.
During multicriteria decision-making, especially when it is applied in the service network design, shippers' preferences are of great importance. Zeng et al. [24] considered the customer preferences and spatial interaction and provided a bidimensional evaluation method to help people understand the influences of Carat Canal. With many studies adopting shippers' preference in freight transportation network optimization, some have noticed the heterogeneity among freight service preferences [25,26]. One study proved that the service level of the network could be improved by considering the shipper heterogeneity [26]. Zeybek [27] highlighted the meaning of designing differentiated segment strategies for rail transport and used a multimethod to distinguish six behavioral customer segments. Liu et al. [28] pointed out the influence of mass customization on the key decisions made by LSPs and integrators in a logistics supply chain. Given that the level of expectations might vary from segment to segment of the customers, models with disaggregated data integrating different factors into decision-making are badly needed.
To the authors' knowledge, few disaggregated models have been applied to the rail transport, especially to the freight network under BRI, even though they are well-researched. Interaction between service network and mass customization could be found because individually designed services require a higher cost [29]. If the LSI could find cooperation with the LSPs, then the service quality would be improved by reducing transportation with low efficiency. Thus, the carbon emission of the network could be reduced. Table 1. Literature using multicriteria decision-making and criteria included.

Price Time Reliability Frequency Flexibility Traceability Emission Safety
Yu et al. [14] García-Menéndez and Feo-Valero [30] Jiang et al. [22] Arencibia et al. [10] Duan et al. [26] Moschovou and Giannopoulos [31] Shinghal and Fowkes [13] Tsamboulas and Moraitis [32] Brooks and Trifts [33] Bergantino and Bolis [34] Reis [35] Zotti and Danielis [36] Wang and Yeo [21] Zeybek [27] Zamparini et al. [37] Norojono and Young [38] Punakivi and Hinkka [39] Wen et al. [40] Zotti and Danielis [36] Beuthe and Bouffioux [41] Witlox and Vandaele [42] Jain et al. [43] Tuzkaya and Önüt [15] Witlox and Vandaele [42] This study aims to explore the decisions regarding achieving an efficient network of CRexpress with mass customization service. When designing the CRexpress service network, the LSI determines an optimal configuration in terms of terminal locations and allocates the shipper demand to different LSPs. This study aims to answer the following three questions: (1) Considering the satisfaction of shippers, what layout plan will the LSI come up with to improve the performance of the LSPs and the overall service network of CRexpress? (2) What type of logistics characteristics should be taken into consideration to capture preferences of shippers? (3) How will the heterogeneous preferences of shippers influence the selection of the consolidation centers, the market share of the logistics services provided by the LSPs, and the decision of the LSI?
The remaining part of this paper is organized as follows: Section 2 describes the system and the problem and then introduces the notation and the bilevel programming formulation. The solution approach to the bilevel model is also introduced. Section 3 includes the application and results. The conclusions and implications can be found in Section 4.

Background of the CRexpress Service Network
The CRexpress service network in our paper is denoted by G (N, A), where N represents the set of nodes (including shippers and gateways), while A stands for the set of arcs in the network. We use the selected 27 cities for the 27 provinces in China in Zhao et al. [23] as the shippers for the CRexpress service network. According to the authors, the cities could be classified into two groups: (i) one group of 17 cities that will not be chosen as consolidation center candidates after evaluation, denoted by I (indexed by i); and (ii) one group of 10 cities as candidates for consolidation centers, represented by J. We further divide the second group into two sets depending on whether it is chosen after optimization. Specifically, city j ∈ J is the optimized consolidation center, and city n ∈ J is excluded from the result. City n ∈ J, which is not chosen as the consolidation center among the candidate cities (set J), will be treated as origin cities in the further analysis. All shippers above have the demand to be exported through the CRexpress network. The railway and road are considered in the CRexpress network, represented by M (indexed by m). Figure 1 is an illustration of the routes and modes of all shippers, which could be summarized as follows: 1.
In every i ∈ I , the cargo is delivered to j ∈ J by railway or road transportation and consolidated into one single CRexpress train with other cargos received at j ∈ J and then continuing to gateway k ∈ K; 2.
In every n ∈ J that is not chosen as the consolidation center, the cargo is delivered to j ∈ J by railway or road transportation and consolidated into one single CRexpress train with other cargo received at j ∈ J and then continuing to gateway k ∈ K; 3.
The transport demand of every j ∈ J will be consolidated with other cargos received at j ∈ J ; then, the cargos continue to gateway k ∈ K; 4.
All cargos are assumed to be transported to the single destination in this study. The shippers' aim from route and mode selection is to identify the optimal alternatives with the highest utility from all feasible alternatives in the network. The utility function is the weighted sum of eight logistics characteristics to evaluate the route and mode alternatives in the CRexpress network. The details of the utilities will be described in the next section.

Utility Function
In this section we introduce the utility function as a foundation for the optimization model. The utility function is the weighted sum of logistic factors to evaluate the alternatives [10]. With the assumptions we build the utility function with integrated indicators in order to improve the decision-making process.

Assumptions
The utility function we build is based on the following assumptions: • All decision makers are economists, who choose the alternative that maximizes their own utility.

•
Only road and railway transportation are included in this system.

•
The decisions made in this system are considered for a specified time period.

•
The indicators in the models are considered constant parameters rather than stochastic variables for the given time period;

•
This system has no congestion, and the loading and unloading times are not considered.

Formulating the Utility Function
We first explain the formulation of city i ∈ I, and the utility function for city n ∈ J could be obtained following the same principle. In each city i ∈ I, the transportation distance consists of two arcs, namely, from the i ∈ I to j ∈ J and from j ∈ J to the destination through k ∈ K. The railway and road are considered, represented by M (indexed by m). The factors cost, time, and emission, which are distanceand mode-dependent, also have two parts. The cost, time, and emission for i ∈ I exporting the cargo equal the sum of the costs incurred by the arcs contained (Equations (1)- (3)).
The factors, such as frequency, reliability, safety, flexibility, and traceability, are only mode-dependent in our case, which means that the second arc has no influence on it. In each city j ∈ J that is selected as the consolidation center, the shippers shall choose from three gateways, and the corresponding factors are listed in Table 2. The measurement scales of the above factors are different, which causes infeasibility to simply add them together. We normalize these attributes for a fair comparison.ĉ m ijk ,t m ijk ,ê m ijk ,r m ijk ,f m ijk ,ŝ m ijk ,p m ijk , andv m ijk are the normalized cost, time, emission, reliability, frequency, safety, flexibility, and traceability of the alternative, respectively. Equation (4) is the min-max normalized function [43], in which we utilize the notation with a line on the top representing the maximum value of the attributes, while the one with a line in the bottom refers to the minimum. All factors are valued in [0, 1] after normalization.
The unit utility function is formulated in Equation (5), where w c , w t , w e , w r , w f , w s , w p , and w v are the non-negative weight parameters for cost, time, emission, reliability, frequency, safety, flexibility, and traceability, respectively. The unit utility function is generally composed of two parts: (i) the negative component of the weighted sum of cost, time, and emission; and (ii) the positive component of the weighted sum of the reliability, frequency, safety, flexibility, and traceability. The weights are the preferences of the shippers toward various attributes when choosing the alternative.
Considering that the total demand could be split into several alternatives, the overall utility of the city i ∈ I is the sum of utility (U m ijk ) multiplied by the volume (x m ijk ) (Equation (6)). The utility of cities n ∈ J and j ∈ J could also be obtained following the same principle (Equations (7) and (8)).

Mathematical Formulation of the Bilevel Model
This model deals with the contradiction among the LSI, LSP, and shippers using the CRexpress network. The leader-follower problem can be mathematically expressed as a bilevel program. The decision-maker in the upper level, the LSI, makes their decision about how to assign the demand flow to different LSPs for minimizing the construction cost. Shippers in the lower level will determine the flows on the network that maximize travel utility.
The decision variables corresponding to the decisions made in the model are listed below. We let binary variable O j represent which city to be chosen as the consolidation center, in which O j = 1 if j ∈ J is selected, and 0 otherwise. We denote x m ijk as the volume transported city i ∈ I to consolidation center j ∈ J to destination through gateway k ∈ K by mode m ∈ M. x m njk is denoted as the volume transported city n ∈ J to consolidation center j ∈ J to destination through gateway k ∈ K by mode m ∈ M, when n ∈ J is not chosen to be the consolidation center. x jk is the ictal demand volume at city j ∈ J that will be transported to the destination through gateway k ∈ K.
Every city in the network has its demand to the destination. In each n ∈ J that is not selected as the consolidation center, it will also serve as a city. Therefore, in each i ∈ I, n ∈ J, and j ∈ J, we let D i , D n , and D j represent the demand of nodes i, n, and j that has to be served by the CRexpress network, respectively. Ca is the maximum operating capacity of a consolidation center. The capacity for each CRexpress train leaving j ∈ J for k ∈ K is denoted by B. c 0 is the fixed construction. The reference of notations can be found in Table A1.
The bilevel problem can be written as follows: Upper level: min Lower level: max m∈M j∈J k∈K m∈M j∈J k∈K m∈M i∈I x m njk ≥ 0, ∀j ∈ J, ∀m ∈ M, ∀k ∈ K, ∀n ∈ J, n j, x jk ≥ 0, ∀k ∈ K, ∀j ∈ J.
In the upper level, the objective function (9) is to minimize the construction cost for setting up a consolidation center. Constraint (10) shows that a consolidation center has to be chosen to be in the network.
In the lower level, the objective function (11) is to maximize the total satisfaction of all shippers when certain consolidation centers are selected by the decision-maker in the upper level (i.e., central LSP). Constraint (12) ensures that all weights add up to one. Constraints (13)- (15) ensure that the demands of all cities are satisfied. Constraints (16)- (17) imply that the flow will not enter any nodes that are not chosen as consolidation centers. Constraint (18) guarantees that the volume received at each consolidation center is under its maximum capacity. Constraint (19) restricts that all cargo leaving each consolidation center are transported by CRexpress. Constraints (20)-(22) define the non-negativity of the decision variables.
In summary, the CRexpress network optimization problem is mathematically expressed as a bilevel program. The upper level objective involves the decision variables in the lower level. We show how to solve this problem in the next section.

Solution Approach and Reformulation
In a given O j from the upper level, the lower level problem (i.e., Equations (11)-(22) is a linear programming formulation. Therefore, the bilevel program we proposed could be reformulated into a single-level MILP formulation by replacing the lower level optimization program with the Karush-Kuhn-Tucker (KKT) conditions [44]. First, we replace the lower level problem with KKT conditions, then, we transform the nonlinearities of the complementarity constraints, finally, we strengthen them to improve the computational performance of solving the program.

KKT Conditions for the Lower Level
To obtain the KKT conditions for the lower level program, first, we describe the primal feasibility together with the dual feasibility. The primal feasibility includes the following, and the dual variables for each constraint in the lower level are parenthetically indicated: m∈M j∈J k∈K x m njk ≤ O j D n , ∀j ∈ J, ∀m ∈ M, ∀n ∈ J, n j, ∀k ∈ K λ m njk , (27) m∈M i∈I k∈K m∈M i∈I Dual feasibility implies the following: Then, we can derive the stationarity equations: The set of complementary conditions is as follows: λ m njk x m njk − O j D n = 0, ∀j ∈ J, ∀k ∈ K, ∀m ∈ M, ∀n ∈ J, n j, ω m njk x m njk = 0, ∀ j ∈ J, ∀k ∈ K, ∀m ∈ M, ∀n ∈ K, n j, µ jk x jk = 0, ∀k ∈ K, ∀ j ∈ J.

Reformulating Complementarity Conditions
The above KKT conditions for the lower level are all linear except for the complementarity set, which can be linearized with a binary variable and a large constant [45]. On the basis of the methods proposed in the literature, we introduce binary variables for constraints (40)- (46), which are the ones with a star on the right top corner (i.e., π m ijk * is the binary variable for π m ijk ) and M as the arbitrarily large number. Then, we obtain the equivalent linear constraints: λ m njk ≤ Mλ m njk * , ∀j ∈ J, ∀k ∈ K, ∀m ∈ M, ∀n ∈ J, n j, m∈M i∈I k∈K m∈M i∈I ω m njk ≤ Mω m njk * , ∀j ∈ J, ∀k ∈ K, ∀m ∈ M, ∀n ∈ J, n j,

Strengthen the Constraints
The large value for M might lead to poor bounds from linear relaxation. In this step, we strengthen some constraints by replacing each big M with a relatively small number, which does not cut off the optimal solution. For instance, in constraint (48), D i is always greater than x m ijk − O j D i ; thus, we can tighten constraint (48) with D i . The rest follows the same principle.
x m njk − O j D n ≤ D n 1 − λ m n jk * , ∀j ∈ J, ∀k ∈ K, ∀m ∈ M, ∀n ∈ J, n j (62) m∈M i∈I k∈K x m njk ≤ D n 1 − ω m njk * , ∀j ∈ J, ∀k ∈ K, ∀m ∈ M, ∀n ∈ J, n j (66) Our proposed bilevel programming model can be finally reformulated as an MILP formulation as follows through the above steps: O, x, y, β, α, π, λ, δ, ρ, η, τ, ω, µ, The final formulation can be accurately solved by existing optimization solvers. We adopt CPLEX 12.9 [? ] to implement this model and perform other scenario analysis, which is discussed in the next section.

Results
In this section, we perform the proposed model in Section ?? by using realistic data as input and analyzing some results for insightful findings. We first describe the data obtained from existing literature. Then, we provide the definitions of all the attributes included and their values. Furthermore, we introduce two problem scenarios to understand the effects of incorporating shippers' heterogeneous preferences into the CRexpress network. Finally, we adopt CPLEX 12.9 as our optimization solver to solve the above model.

Input Parameters
In this case, we utilize the network as used in Zhao et al. [? ]. Twenty-seven cities are present in the network, among which 10 are chosen as alternative consolidation centers. Considering all the three possible corridor routes, we have Manzhouli, Ernhot, and Altwa as gateways. We assume Hamburg to be the only destination. The annual volume between each city and Hamburg in 2020 is calculated on the basis of Zhao et al. [? ]. The maximum capacities of the consolidation center and CRexpress trains are 10,000 and 2200 t, respectively. The fixed cost for setting up a consolidation center is 6.8 × 10 8 .. Now, we describe the data for all the logistic attributes, namely, price, time, frequency, reliability, flexibility, traceability, safety, and emission.
According to the China Railway Customer Service Center website (https://www.12306.cn), three types of services are offered, namely, full truck load, less than truck load, and container transportation. We follow the method in Zhao et al. [? ] to obtain the average price. First, we define the maximum container load of the common international standards for containers as 22 and 27 t for 20-and 40-foot containers, respectively. On this basis, the average freight train cargo equals The final formulation can be accurately solved by existing optimization solvers. We adopt CPLEX 12.9 [46] to implement this model and perform other scenario analysis, which is discussed in the next section.

Results
In this section, we perform the proposed model in Section 2 by using realistic data as input and analyzing some results for insightful findings. We first describe the data obtained from existing literature. Then, we provide the definitions of all the attributes included and their values. Furthermore, we introduce two problem scenarios to understand the effects of incorporating shippers' heterogeneous preferences into the CRexpress network. Finally, we adopt CPLEX 12.9 as our optimization solver to solve the above model.

Input Parameters
In this case, we utilize the network as used in Zhao et al. [23]. Twenty-seven cities are present in the network, among which 10 are chosen as alternative consolidation centers. Considering all the three possible corridor routes, we have Manzhouli, Ernhot, and Altwa as gateways. We assume Hamburg to be the only destination. The annual volume between each city and Hamburg in 2020 is calculated on the basis of Zhao et al. [23]. The maximum capacities of the consolidation center and CRexpress trains are 10,000 and 2200 t, respectively. The fixed cost for setting up a consolidation center is 6.8 × 10 8 . Now, we describe the data for all the logistic attributes, namely, price, time, frequency, reliability, flexibility, traceability, safety, and emission.
According to the China Railway Customer Service Center website (https://www.12306.cn), three types of services are offered, namely, full truck load, less than truck load, and container transportation. We follow the method in Zhao et al. [23] to obtain the average price. First, we define the maximum container load of the common international standards for containers as 22 and 27 t for 20-and 40-foot containers, respectively. On this basis, the average freight train cargo equals 0.15 RMB/ton/km. With regard to the road freight rate, we use 0.45 RMB/ton/km, which is calculated as the price estimated from a Chinese road freight website (https://www.jctrans.com).
The traffic safety regulations in China set the maximum and minimum speeds for trucks on highways to be 100 and 60 km/h, respectively, considering different situations in China, the speed of 60 km/h is applied for road transportation. The average speed of freight trains in China has been estimated to be 40 km/h, and the speed of CRexpress is defined as 60 km/h because the speed within China is slower than those overseas. Distance data are obtained from the China Railway Customer Service Center website (http://www.12306.cn) and the China highway website (http://www.china-highway.com). The transportation time for road and railway inside China could then be calculated (Tables A4-A6).
We use the average frequency of two trains/day for railway inside China. In real life, road trucks can leave at any time. However, a large value may introduce numerical stability issues, thus, we suppose that road trucks leave every hour. Based on the data released on https://www.yidaiyilu.gov.cn/xwzx/ gnxw/76434.htm, 6300 CRexpress trains operated in 2018. Considering the diversity in different parts of China and three corridors, we utilize three trains per day as the average frequency of CRexpress.
The attribute reliability is regarded as the percentage of cargos transported to the destination within a given time [35]. Safety refers to the percentage of the commercial value of the total freight shipped that is not lost or damaged; while transportation flexibility represents the percentage of executed nonprogrammed shipment without undue delay [41]. Railway is the more reliable than road transportation [35]. For example, road trucks can be easily disturbed by traffic conditions or bad weather, while freight trains are less affected. With regard to transportation security, railway could provide more security than road transport. By contrast, road transportation has a good flexibility performance in responding to unexpected demands because rigid timetables are not that necessary and the ease of freight loading and unloading in road transportation [15].
These last years, transportation traceability has become the focus of interest for many researchers and decision-makers because it contributes to the control of flows and allows the optimization of the supply chain in its totality. Many investigations on the system have been conducted to support a real-time access to information for supply chain stakeholders [47,48]. Different shippers may vary in the satisfaction toward the traceability they received from different supply chain suppliers. Therefore, we use percentage in our study to represent shipper satisfaction toward traceability. The Global Positioning System (GPS) and Radio Frequency Identification (RFID) can help provide the shippers additional information about their cargos, however, it comes at a compromise cost. The level of application of GPS and RFID and tracking of cargos is higher in road than in railway/CRexpress. Therefore, shippers will have higher satisfaction toward road than railway/CRexpress. Safety, flexibility, and traceability for road and railway transportation are randomly generated corresponding to the discussion above. Further details about the conversion into percentage can be found in [40].
Greenhouse gas emission is the main environmental criterion for the CRexpress' performance. Carbon dioxide, methane, nitrous oxide, and ozone are the primary transportation-related, man-made greenhouse gases. The impact of these gases can be calculated together as carbon dioxide equivalents (CO 2 e). The model in this study focuses on the comparison between road and railway, and all parameters can be normalized into values in [0, 1]. Accordingly, we use the emission of carbon dioxide (CO 2 ) to represent the emission criteria.
Similar to the method proposed in Wiegmans and Janic [49], we calculate the environmental performance of road and railway transportation with data collected from Jong and Riet [50]. Road transportation generates a higher level of carbon dioxide than freight trains. Table 3 summarizes the values of eight criteria for road and railway inside China.

Scenarios Description
In the application of the proposed bilevel model, two scenarios are included in which we have two representatives of shipper population. In scenario A, all shippers hold homogeneous weights, whereas in scenario B, shippers distinguish with heterogeneous preferences. Unlike in scenario A, where all shippers share the same preference, each demand flow in scenario B is divided into three clusters from that in scenario A on the basis of the membership size. The different preferences and sizes of three clusters are obtained from Li et al. [51].In the aforementioned study, a stated preference survey was conducted from 29th March 2019 to 18th April 2019 via face-to-face personal interviews at Chengdu International Railway Service Co., Ltd. in China. This company is one of the important platforms where CRexpress' customers from all over the country gather. Therefore, the result of their study could represent the preferences of all CRexpress' customers in China. Table 4 lists the different preferences and cluster sizes. The visualized weighted, directed graph is shown in Figure 2. Scenario A simulates the baseline where we have the aggregated weights of the attributes. This scenario could be concluded from that in the baseline scenario, and the weights for all criteria are relatively equal, which is quite different from that in scenario B. Cluster 1 (38.1% of shippers) consists of shippers who are highly sensitive to price, followed by time and reliability, and the weights of rest criteria are rather small. The shippers in cluster 2 (50.8%) are also easily affected by price, however, reliability also means much to them. With regard to the shippers in cluster 3 (11.1%), the important indicator is reliability, and their preferences for price are not as strong as the former two clusters.
If we compare the weights of three clusters in scenario B with the aggregated weights in scenario A, cluster 2 is the most similar one, while clusters 3 and 1 are quite different. This finding shows that the traditional way of modeling ignored the type of shippers in clusters 1 and 2, which are observed in our study.

Designs and Analysis
We obtain the optimal layout of flows and the selection of construction centers in two scenarios, namely, designs (A) and (B), by applying the proposed model to the above inputs. Taiyuan, Zhengzhou, Wuhan, Chengdu, and Suzhou are chosen in design A, while design B identifies Xian, Lanzhou, Taiyuan, Wuhan, and Suzhou. The details of the flows, including the transportation mode from cities to consolidation centers in two designs, can be found in Table A2. Details of flow are in design A and Table A3. Most cities generally change their choices either for transportation mode or consolidation centers. We then investigate the adjustment from two aspects, namely, (1) the whole shippers and (2) difference among three clusters.
The comparison from the perspective of the whole network is based on the fair calculation by using the heterogeneous preference in scenario B, given that the three types of shippers exist in real life. Table 5 illustrates the utilities of all cites in both designs. The comparison result of the utilities between the two designs indicated that 15 out of 27 cites increased in utilities (indicated with up arrows) after segmentation, with a volume of 31,991 t, which is 78% of the total demand volume. Cities with an underline are the potential consolidation centers. Given that the utility of cost, time, and emission are negative, the total utility of some cities might be negative (Table 5). ↑ and ↓ represent the increase or decrease of utilities, respectively.
To achieve a highly visualized impression of the key performances in two designs, we introduce the term of performance for the eight indicators involved. The performance of each indicator is measured with the sum of flow volume multiplied by its corresponding utility. The total performance is the sum of the eight indicators. Figure 3 presents all the performances, including the total performance and eight indicators. Given that the weight of emission, time, and cost are negative, their unit performances are also negative. The negative indicators that are closer to zero represent the better performance level. Meanwhile, the positive ones with higher values are highly preferred.
The total performance level is improved because the total performance in design B is higher than that of design A by 10.4%. More than half of the shippers have their service quality improved, and the number of consolidation centers are the same in two designs, hence, including shippers' heterogeneous preferences will not bring any extra construction cost. Accordingly, the model with shippers' heterogeneous preferences greatly improves in providing services that could satisfy many shippers. Furthermore, the analysis in Figure 4 shows that the share for railway as the transportation mode from origin city to consolidation center increased by 13% in design B than in design A. Railway transportation has the advantage of low emission, which is only 20% of road transportation. Therefore, the model with segmentation not only improved the service level but also provided a sustainable network. Different clusters have varying reactions to the segmentation.  Figure 5 demonstrates the utilities of three clusters in two designs, with utilities of all shippers on the right. The shippers in cluster 1 have huge promotion, those in cluster 3 have decreased utility, and those in cluster 2 suffer most. This phenomenon occurs because in design A, we ignored the type of shippers in cluster 1 ( Figure 2). However, the utilities of such shippers are increased by recognizing them in design B. The increase is the compensation made by the shippers in clusters 2 and 3, who share similar preferences in terms of aggregate weights. The shippers in cluster 3 have different priorities in reliability and price with the aggregate weights, which is shared by shippers in cluster 2. Accordingly, sacrifice of these shippers is slighter than those in cluster 2. In conclusion, Figure 5 confirms that the improvement of utilities of all shippers in design B is achieved by recognizing the different types of shippers (cluster 1), which is one of the contributions of our work.   Figure 5 illustrate that shippers with less utility after segmentation are mainly the underlined cities, which are defined as potential consolidation centers in set J. Although some cities, such as Shijiazhuang and Hefei, which are not in set J, also suffered, their sacrifice are less than those in set J such as Taiyuan and Chengdu. This notion means that the decrease observed in Table 5 is the result of the decrease in set J.
The details in Table A2. Details of flow in design A provide the explanation of the above findings. The cities in the potential consolidation centers (set J) choose farther gateways if they are chosen. If a city in set J is not chosen as a consolidation center, then the demand is likely to be transported in the most satisfying way for the shippers. Nevertheless, design B provides the priority of the nearest gateway to the volume from other origin cities if the city is chosen as a consolidation center because of the capacity of CRexpress. The chosen consolidation centers must satisfy the demand of unchosen origin cities to improve the overall level of the whole system.
The following analyses are conducted on two groups, i.e., cities in sets I and J, for fair calculation. Two reasons have to be stressed for group analysis. First, city n in set J will be treated as an origin city if it is not chosen as a consolidation center. If we divide the cities on the basis of whether it is an origin city, then the total volume will change. Second, most cities suffered if they are chosen as consolidation centers, which means that the choices are a combination of their real preferences and the purpose of optimizing the network. Accordingly, the impact of shipper real weights on the choice is difficult to investigate.
The analysis of mode share of shippers in set I confirms that the different reactions among three clusters are in line with their importance toward each logistics character in Table 4.
We first calculated the mode share of the three clusters in two designs ( Figure 6). Shippers in clusters 1 and 2 show more preference for railway in design B than in design A, while those in cluster 3 have the opposite choice. Table 3 illustrates that railway has the advantage of price and emission compared with road. Accordingly, shippers who show strong importance toward price and reliability are likely to choose railway. However, the weight for price in cluster 1 increased by 58% compared with that in design A, while the rate for cluster 2 is only 38%. The weight for reliability increased by 58% for shippers in cluster 3, while that for price decreased by 41%. Based on the discussion above, the different mode shares of shippers in set I are quite reasonable. The overall share of railway increased by 13%, and shippers with high weight for price and reliability can get their cargo transported by railway in design B. Therefore, our proposed model makes a great contribution to improving the sustainability of the system because it meets most shippers' demands. We summarize the eight indicators among three clusters in two groups of cities. The results of indicators can be divided into two parts, namely, the positive (Figure 7a,b) and the negative (Figure 7c,d). We summarize all the negative indicators among three clusters, also with that from design A as the baseline to the left. The heights of each indicator could tell us how much it has been improved before and after segmentation. The eight indicators are compared to see how the preferences of shippers influence the decision-making process, thus, we focus on the factors that are expected to be improved in design B. Specifically, shippers in cluster 1 expect a lower price, those in cluster 2 desire lower price and high reliability performance, and those in cluster 3 would prefer a high reliability and safety level ( Figure 2).
First, we look at the cities in set I. Apparently shippers in clusters 1 and 2 get exactly what they want, while those in cluster 3 do not. This situation may happen because the size of cluster 3 is only 11.11% of all the shippers in set I, which could be easily made to compensate for optimization of the whole system. With regard to cities in set J, shippers in cluster 1 received the lowest price, just as their preferences. Shippers in cluster 2 still have high reliability, but it comes with a high price. This phenomenon occurs because such shippers sacrifice the further gateway for convenience of other origin cities. Shippers in cluster 3 have increased reliability and safety, however, regardless of which gateway the chosen consolidation centers choose, the normalized values of reliability and safety are always one, which is larger than that of the origin cities, which vary in [0, 1]. The increase in cluster 3 does not mean that shippers are more satisfied in design B.
The model built after segmentation on the basis of preference heterogeneity could provide a service network for CRexpress with high service quality and sustainability level for most shippers in three ways. First, although some shippers have decreased utility, the utility of the total system can still be increased, which is exactly the purpose of the LSI. Second, the total share of railway transportation, which is a more sustainable mode compared with road, is promoted by 13%. Finally, the LSI of the CRexpress network could provide the shippers with additional customized services with the integrated indicators, which will make a great contribution to the maintenance of its competitiveness.

Discussion and Conclusions
The service network design problem with customization is an important problem, which deals with the contradiction among the LSI, LSPs, and shippers. The traditional way of operating the network classifies the cargo on the basis of its size and appearance, which does not pay great attention to shippers' preferences toward other logistics characteristics, and is no longer suitable at this stage. The establishment of the BRI together with the increasing international trade between China and European countries also places an additional requirement for the CRexpress to maintain its competitiveness. A bilevel programming formulation is established to address the necessity of integrating logistic service-related attributes and shippers' heterogeneity in the CRexpress network.
In the upper level of the model, the construction cost is minimized to reduce unnecessary transportation, while the lower level maximizes the satisfaction of shippers. The model is transferred through the methodology of the KKT conditions into a linear programming formulation and then solved with an advanced commercial optimization solver. We identify eight indicators, namely, price, time, reliability, frequency, safety, flexibility, traceability, and emission, to assess the performance of the network by reviewing the literature with integrated attributes.
Two designs are established in the case study to demonstrate the impact of considering preference heterogeneity, i.e., (A) based on preference homogeneity and (B) divided into three groups with different preferences. The comprehensive analysis between two designs indicated that the utility of more than 78% of the shippers improved while maintaining the same construction cost by considering the variations between shippers in the optimization process. Many shippers choose railway as the transportation mode as a result of segmentation, which is also the desire of the LSI because it is sustainable. Further comparison of eight indicators shows that most shippers could get the most promotion in the criteria that they weighted most highly. The above results imply that shippers' heterogeneity may be allowed to improve shippers' satisfaction to conduct a well-customized service network.
This study could also serve as a reference on how to provide proper subsidies to the cities' operation of CRexpress. The knowledge about the preferences of different user groups and the performance of different designs can further be used by the government to establish how subsidies should be allocated to achieve the desired effects. For example, if the objective of a government subsidy scheme is to achieve a greener transport system, the emission reduction obtained with a design, combined with an assumption about market prices, can be used to evaluate the required subsidy to achieve a unit value of emission reduction in the network. Subsidies could be allocated based on these costs.
In summary, this study contributes to the optimization of the CRexpress network with four aspects. First, the bilevel model builds a good connection among the LSI, LSPs, and shippers. This aspect not only gives the selection of consolidation centers but also the optimal assignment of the flows of all shippers to different LSPs. Second, the proposed model is reformulated by the KKT conditions and solved with the commercial optimization solver, which is of great interest to operators in real life because it is simple and easily accessible. Third, the model contains eight different logistics characteristics, which could provide a comprehensive description of the system. Lastly, shippers' heterogeneity is effectively considered in the model and elicits significant improvement.
The potential directions of the further study of this topic include but are not limited to the following topics. First, further study should be extended to the impact on CRexpress of other transportation modes, such as shipping and airfreight, given that China has become increasingly important in international freight transportation. Second, a game theory model should be taken into consideration to explicitly understand the cooperation and competition among cities, together with that between the shippers and the LSP.

Conflicts of Interest:
The authors declare no conflict of interest. Table A1. Notations and definitions.

M
Set of transport modes, indexed by m I Set of cities, indexed by i J Set of cities that could be chosen as consolidation centers, indexed by j if it is chosen, and n otherwise K Set of gateways, indexed by k Parameters D i Demand volume from city i D n Demand volume from city n D j Demand volume from city j U m ijk Route utility from city i to consolidation center j to the destination through gateway k by mode m U i , U n , U j Total utility of all shippers in city i, n, and j, respectively c o Fixed construction cost for setting up a consolidation center B Maximum capacity of a single CRexpress train Ca Maximum operation capacity at a consolidation center k M A large-enough number set by the user Decision variables x m ijk Volume transported from city i to consolidation center j to destination through gateway k by mode m x m njk Volume transported from city n to consolidation center j to destination through gateway k by mode m x jk Original demand volume at consolidation center j to the destination through gateway k O j O j = 1 if consolidation center j is selected, and 0 otherwise