#### Abstract

Congestion pricing is one effective demand management strategy to alleviate traffic congestion. This work investigates pricing schemes for mixed traffic flow systems where the human-driven vehicles (HVs) and autonomous vehicles (AVs) coexist. The emerging and integration of autonomous vehicles can help improve the overall transportation efficiency and safety. Given the coexistence of HVs and AVs in the near future, there is need to adjust the existing traffic management strategies to adapt to the mixed traffic conditions. In this study, congestion pricing is imposed on the HVs and the AVs differently, that is, a distance-based toll to the HVs while a delay-based toll to the AVs. We consider six user groups based on the value of time (VOT) and the vehicle types. Compared with the unified distance-based toll, the advantage of delay-based toll is demonstrated first. Then, a surrogate-based optimization framework, namely the regressing Kriging (RK) model, is formulated. Three pricing schemes are investigated and compared: equity-oriented (EQ), environment friendliness-oriented (EN), and revenue-oriented (RE) schemes. Results show that the RE scheme collects the highest revenues; however, its cost-efficiency is weakened. The EQ scheme reduces the variance in the average travel costs among user groups, thus solving the equity issue.

#### 1. Introduction

With the rapidly increasing population, traffic congestion has become a major problem, especially in urban areas. Due to the limited land resources and high costs, it is unsustainable to address this issue by merely building new infrastructures. Given this situation, congestion pricing, one of the promising demand-oriented strategies, is widely advocated to mitigate congestion. As an economic lever, it originates from the Pigou’s theory where extra traveler entering the network should be charged an additional fee for his negative impact imposed on other road users [1]. This type of pricing model is known as the first-best pricing. However, the first-best pricing is difficult to be implemented because it is unrealistic to charge all links in a network strictly by the marginal cost surplus of traffic [2]. Instead, the second-best pricing model was proposed wherein just a subset of links are charged based on the traffic conditions [3]. A branch of congestion pricing practices emerged since 1975, when the area licensing scheme (ALS) was successfully operated in Singapore. However, quite a few schemes were aborted at the stage of trial or referendum. Low public acceptance was the main reason which led to the failure of these cases [4].

Equity is a vital factor influencing public acceptance towards the congestion pricing policy [5]. As a result, vast literature discussed the equity issue on the design of congestion pricing schemes from different perspectives. Lucas et al. pointed out that drivers’ places of residence will affect their willingness to pay. Under a cordon-based pricing scheme, people living outside the cordon and commuting towards it may become losers while the winners are those living and conducting activities inside the cordon [6]. Gaunt et al. investigated residents’ attitudes towards congestion pricing in Edinburgh (Scotland). He found that people with high car dependency were less supportive of congestion tolls [7]. Chen at al. conceived a novel pricing scheme considering drivers’ perceived level of service (LOS). People who enjoyed higher LOS on their trips would pay more [8].

Value of time (VOT) is an important direction when we deal with the equity issue. According to Arnott’s research, the heterogeneity in VOT affected the calculation of welfare effects. The congestion toll favored people with high VOT relative to the schedule delay [9]. Later, more studies were dedicated to examine the impact of congestion pricing on users with different VOTs [10–12]. However, few researches designed the equitable toll rates directly based on the difference in VOTs. Zheng et al. carried out a seminal study where the group-based tolls were determined according to the distribution of VOTs. They argued that when the total time savings after the implementation of congestion pricing were monetized, drivers with higher VOT would gain more. Hence, it is unfair to ask all road users to pay the same amount of toll [13]. However, the VOT categorization in their research is only related to income levels. Obviously, the vehicle type also influences VOT. The advent of autonomous vehicles (AV) has the potential to impact travel behaviour as well as the transportation system. The AV technologies enable drivers to do other activities on their trips, thus making travel time more productive [14–16]. Zhong et al. applied a mixed logit model to quantify the changes in VOT if taking AVs. The results showed that the VOT would be reduced by 30% at most [17]. On the contrary, the majority of actualized congestion pricing schemes are cordon-based or distance-based [5]. The drawback of the cordon-based scheme is that drivers are charged equally regardless of the actual distance travelled inside the cordon. As for the distance-based case, it may induce people to choose the shortest paths, thus making the distribution of congestion more uneven. Gu et al. once put forward the joint distance and delay toll (JDDT) scheme to remedy the limitation of the distance-based congestion pricing policy [2]. But the novel scheme is difficult to be put into practice in the era of human-driven cars. With the embedded vehicle-to-vehicle communication technology, AVs can collect rich information for the real-time traffic condition. As a result, it provides opportunities to apply advanced congestion pricing approaches.

The traditional way to solve the congestion pricing problem is formulating a bi-level mathematical program with equilibrium constraints (MPECs) [18–20]. Nonetheless, it remains challenging to deal with MPEC in a dynamic large-scale transportation network. As the dynamic congestion pricing (DCP) problem often involves high-dimension decision variables and a nonconvex, nonlinear, and nonclosed form objective function [21–23]. In recent years, with the help of mature dynamic traffic assignment (DTA) simulators, simulation-based optimization (SBO) is recognized as an alternative way to handle this problem. It does not require an explicit mathematical formulation of the objective function which is usually unavailable under stochastic traffic dynamics. In general, SBO methods can be classified into four types: (1) direct search; (2) gradient-based approach; (3) feedback control; (4) surrogate-based method. The direct search partitions the search space into hyper rectangles. By comparing the objective function values, the potential optimal rectangle is determined for further partition [24]. Given its exhaustive partitioning characteristics, it is time-consuming to improve the current best solution. The gradient-based method is derivative-based, which tries to approximate the gradient by finite difference. However, it is not suitable for objective functions where the gradient does not exist everywhere [25]. The feedback control is an efficient optimization method. It adjusts the input variables iteratively to let the system output approach the set point [26]. But this method is not applicable when it comes to high-dimension problem with complex constraints [27]. In this context, the surrogate-based methods gain popularity for their ability in approximating expensive-to-evaluate functions. With limited number of evaluations, they can approximate simulation input-output mapping [28]. The regressing Kriging (RK) model outperforms other surrogate models for its excellent performance in prediction [29]. By constructing a probability model, it generates not only an interpolated spatial correlation but also an estimate of the uncertainty. Meanwhile, the feature of regressing sample points overcomes the influence of simulation noise [30].

In order to monitor the effects of the congestion pricing schemes, it is essential to understand the traffic dynamics at the network level. Geroliminis and Levinson proposed the macroscopic fundamental diagram (MFD), which related the network average density to the average flow [31]. The congestion can be identified when the flow decreases with the density, and it indicates the oversaturated state of the network. The critical density is the one at which the maximal production (flow) is achieved [32]. After the implementation of the pricing scheme, we hope the network density would be around the critical value.

Although plenty of studies attempted to design appropriate congestion pricing schemes, there are still some research gaps: (1) few works considered the mixed traffic of human-driven vehicle (HV) and AVs. The inserted communication technologies of AVs enable us to apply the advanced pricing strategy, which can be different from that of HVs; (2) reduction in VOT for AV users was not well considered; and (3) existing studies rarely compared the pricing schemes with different objectives. In this paper, we propose a surrogate-based optimization to solve the dynamic pricing problem for mixed traffic with drivers of heterogeneous VOTs. The HVs are charged by distance while AV users pay according to the average link delay. The VOT-based tolls are adopted among different groups. The results of three pricing schemes, namely equity-oriented (EQ), environment friendliness-oriented (EN), and revenue-oriented (RE) schemes are compared. To the best of the authors’ knowledge, the congestion pricing optimization under mix traffic with heterogeneous users has never been proposed before and hence offers significant contributions to the literature.

#### 2. Methodology

##### 2.1. Simulation-Based Dynamic Assignment

The stochastic route choice (SRC) model is a kind of nonequilibrium-seeking DTA model. It has been extensively used for transportation operations owing to the lower computational cost compared to the dynamic user equilibrium (DUE) model [33]. Table 1 summarizes the notations used in SRC model. Let denote the object network. The time-varying distance-based toll rates for the HV user group can be represented as a decision vector . And the delay-based toll rates for AV user group is expressed as .

The generalized travel cost for HV user group to choose link in the interval is as follows:

The generalized travel cost for AV user group to choose link in the interval is as follows:

Note that the communication is only among AVs. Hence is the average value acquired from AVs which crossed link during the time interval. Given the mixed traffic, it can approximately represent the average delay for all groups. The C-logit model is used to simulate the route choice behavior. The probability of choosing path among the path set for a given OD pair in the interval is as follows:where the term () is the commonality factor of path which describes the degree of overlapping with other alternative paths [34]. is the sum of the travel costs of all links in path .

##### 2.2. Proposed Framework

In this paper, travelers are categorized according to the income level and the vehicle type (HV or AV). Given the same vehicle type, high-income travelers have higher VOT due to high productivity [35]. Given the same income level, AV users will have a smaller VOT because they can spend in-vehicle time engaged in other activities. The differentiated VOTs influence the route choice behavior. When paying the same amount of toll, travelers of high VOT will perceive fewer costs (according to equation (1) and (2)), thus resulting in inequity. As a result, let denote the distance-based toll rate for the HV user group of the lowest VOT (HV reference group) in the time interval. Then, the toll rates for other HV user groups are set based on the VOT difference from the HV reference group. Similarly, denote the delay-based toll rate for the AV user group of the lowest VOT (AV reference group) in the time interval as . Then, the toll rates for other HV user groups are set based on the VOT difference from the AV reference group. We focus on the optimization of congestion pricing for the congested central area of the city. Three different congestion pricing schemes are investigated: (1) EQ scheme; (2) EN scheme; and (3) RE scheme.

As we adopt the distance-based and delay-based tolls for HV and AV user groups separately, the average travel costs (AUD/trip) for different groups would be disperse after the congestion pricing. Hence, the objective of the EQ scheme is to minimize the variance in average travel costs:where is the objective function for the EQ scheme and is the average travel cost for the user group after the implementation of congestion pricing. Note that includes the trip travel time (in monetary unit) and toll paid.

As for the EN scheme, we aim to minimize the total emission in the pricing zone. The QUARTET pollution emission model is utilized [36]. In the mesoscopic simulation, we can only acquire the average speed for a link during each time interval. The impact of the acceleration and deceleration processes has been considered via the average speed:where is the objective function of the EN scheme; is the total emission of the user group during the toll period; and is the number of user groups.

The RE scheme is devoted to maximizing the total toll paid. It is equivalent to minimize the negative version of this term. The collected revenues can be redistributed to road users or used to maintain the transport infrastructures:where is the objective function of the RE scheme; is the revenues collected from the user group ; and is the number of user groups.

There are several constraints considered for all three schemes:

Constraint 7 controls the network performance. is the network average density in the time interval; is the critical density of the network; and is the preset threshold. This constraint makes sure that the network density is close to the critical value over time. Constraint 8 and 9 ensure that the toll rates between adjacent time intervals do not change sharply. is the distance-based toll rate for the HV reference group in the time interval; is the delay-based toll rate for the AV reference group in the time interval; and and are the smoothing parameters. Constraint 10 and 11 set the upper and lower bound for the toll rates. is the lower and upper bound for the distance-based toll rate of the HV reference group. is the lower and upper bound for the delay-based toll rate of the AV reference group. Note that the toll rates of other groups are set proportionally. Once the toll rate ranges for reference groups are confirmed, the ranges for other groups are also determined.

##### 2.3. Surrogate-Based Optimization

The RK model is employed to do the constrained optimization. Figure 1 shows the procedure of the surrogate optimization. Some initial points (or toll plans) are generated by a certain sampling method. In our work, Latin hypercube sampling (LHS) is adopted. It is a space-filling approach which stratifies each dimension of the decision variables into an equal number of partitions [37]. At least initial points are needed for -dimension problems [38]. Considering the high dimension of our optimization problems, totally points are generated. The we run the network simulation to get the objective function values for these initial inputs. The simulation input-output mapping is used to build a preliminary surrogate model. In order to enhance the model, additional infill points are required. The additional samples can provide more information in the potential regions where the good solution may exists. The infill process will not terminate until the stop criteria are met. Finally, the accuracy of the model is checked. When the accuracy is validated, we can declare that the optimum is achieved. Otherwise, the surrogate model should be redesigned.

###### 2.3.1. Model Construction

For the RK model, the prediction for a nonsampled point is written as follows:where is the constant mean; is a normally distributed and independent estimation error, and it can be expanded as follows:where is the output vector for sampled points ; is the identity matrix; is the correlation matrix between the sampled points; and is the regression constant. With the positive , the RK model regresses the data to reduce the influence from the simulation noise. is the correlation vector between sampled points and the new point. The correlation between point and is given as follows:where is the dimension of variables; and is the scaling and smoothness coefficient for the dimension, separately; and denotes the sensitivity to the objective function value. A larger means the certain dimension is more sensitive than others. The use of is usually advocated. However, proves to be suitable for the engineering-based problem [39].

Then, the vector of scaling coefficients , the variance and the mean are estimated by maximizing the likelihood probability [40]:

The limited memory Broyden–Fletcher–Goldfarb–Shanna (LBFGS) algorithm is employed to solve the maximization problem. This iterative method is appropriate for problems with large numbers of variables [41].

###### 2.3.2. Infill Strategy

In order to improve the accuracy of the surrogate model by the augmented data set, the acquisition function is used to guide the search direction. The expected improvement (EI) function is a widely-used acquisition function, which estimates the magnitude of improvement at a new point [42]. It can explore the unvisited region and exploit the domain of interest at the same time. Let denote the current best solution. The improvement of a new point can be expressed as . Then, the formulation of EI is demonstrated as [43]:where is the normal cumulative distribution function; is the probability density function; and is the re-interpolation prediction error for the objective:where is the estimate of for the re-interpolation:

Note that by replacing with , stays at zero for all the existing points. Hence, the EI sampling will not be trapped at existing sample points. The new points are obtained with the aim of maximizing the .

The EI function is sufficient for unconstrained optimization. However, when it comes to the constrained optimization. A new point that has a good objective value may not meet the expensive-to-evaluation constraint at the same time. As such, another acquisition, namely the probability of improvement (PI) function is utilized. The mechanism of PI is to research for a new point which has the highest probability to meet the complex constraint. We need to construct a new surrogate model for the constraint separately. Denote as the estimated constraint value. Then, we can formulate the PI function as follows:where is the preset threshold as in Constraint 7 and is the re-interpolation prediction error for the constraint. The constrained EI function can be written as follows:

###### 2.3.3. Model Validation

The leave-one-out cross-validation is a method frequently used to measure the accuracy of surrogate models [44]. This approach leaves out one observation and reconstructs the surrogate model with the remaining points. The estimated output of the removed point from the re-fitted model will be compared with the actual output. The key metric for the leave-one-out cross-validation is the standardized cross-validated residual (SCVR). It is calculated as follows:where is the actual output for ; is the estimated output of without ; and is the prediction error of from the re-fitted model without . The model is validated when the SCVRs for all points roughly fall in the interval [−3,3]. It means the model proves to be 99.7% confident that lies in the interval [, ] [41].

#### 3. Results and Discussion

##### 3.1. Experiment Setup

In this paper, we perform a mesoscopic simulation-based DTA in the Melbourne city, Australia in the traffic simulation software AIMSUN. The travel demand is collected via the loop detector and changes every 15 min between 6 : 00 AM and 10 : 00 AM. Figure 2(a) shows the congested central business district (CBD) (inside the red box). We choose it as the pricing zone. There are 75 nodes and 218 links in the object network. The information of daily income of residents is extracted from the Victorian Integrated Survey of Travel and Activity (VISTA). We assume that people work 8 hours a day. As a result, the VOTs equate to the daily income divided by eight. Travelers are classified into three groups based on the distribution of VOT initially. The VOT of each group is the median of the corresponding population. Since the mass market penetration of AVs is challenging in the near future, 30% AV penetration rate is used in our work. Meanwhile, we assume that more high-income travelers can afford AVs. Finally, there are totally six user groups shown in Table 2. Note that AV users will have a reduction in their VOTs by 30% [14].

**(a)**

**(b)**

**(c)**

The agents of different user groups are uniformly distributed in the network. In order to avoid simulation randomness, we run three replications with different random seeds. The MFDs for the three replications (R1, R2, and R3) and the average result (AVE) are presented in Figure 2(b). According to the MFDs, the flow starts to drop at 25 veh/km/lane. Hence, we set the critical density as 25 veh/km/lane. The network density over time in the no-toll scenario is demonstrated in Figure 2(c). Observe that the density exceeds the critical value around 8:30. So the toll period is from 8:30 to 10:00.

In our study, is set as 0 AUD/km, and is set as 1 AUD/km. is set as 0 AUD/h, and is set as 15 AUD/h. The smoothing parameter is set as , while is set as . One can refer to the study by Gu et al. [34], wherein the analysis on the smoothing parameters is conducted. The threshold is set as 2.45 veh/km/lane which is nearly one-fourth of the value (9.8 veh/km/lane) in the no-toll scenario. The toll period is divided into six 15-minute intervals. Hence, the decision variable vector has 12 dimensions, six of which are distance-based toll rates for the HV reference group, while the rest are delay-based toll rates for the AV reference group.

##### 3.2. Comparison between Unified Delay-Based Toll and Distance-Based Toll

In order to unveil the advantage of the delay-based toll, we compare the network performance under the unified delay-based toll and distance-based toll. First of all, let us introduce a term, namely the spatial spread of density [45]. It measures the heterogeneity of congestion distribution and is calculated as follows:where is the length of link ; is the number of lanes of link ; is the density of link ; and is the network density.

The relationship between the spatial spread of density and the network density can be expressed as the third-polynomial function [34]:

We collect 1000 data points to plot Figure 3(a). The lower envelop of Figure 3(a) is fitted to Equation (19). Parameters are estimated as: , , and . According to Figure 3(a), the spatial spread of density will naturally increase with the network density. Simoni et al. further proposed an advanced metric: the deviation from spread [46]:

**(a)**

**(b)**

Smaller indicates the lower heterogeneity of congestion distribution, thus leading to the higher network flow. We select a delay-based toll scheme with the unified delay toll rate of 13 AUD/h and a distance-based toll scheme with the unified distance toll rate of 0.5 AUD/km. According to Figure 3(b), the network density under the distance-based toll scheme (DI_density) is smaller than that of the delay-based toll (DE_density) in most of the time. However, the result of the network flow is the opposite. When we look into the average deviation from spread during the toll period, we find that it is much smaller under the delay-based toll scheme (1.93 veh/km/lane vs. 4.93 veh/km/lane). This conclusion is in line with Gu’s research [27]. It indicates that the delay-based toll can reduce the heterogeneity of congestion distribution to a great extent. This is also one motivation of our research to implement the delay-based toll on AVs which are equipped with the perfect technology to realize it. The inserted advanced communication technology of autonomous vehicles enables us to collect data such as the vehicle position and the average delay of queues.

##### 3.3. Results of Optimal Toll Designs with Different Objectives

The stop criterion in this paper is that 30 infill points are collected or for consecutive 10 points. Figures 4(a)–4(f) show the results for the EQ scheme. When we track the convergence history of EI, the value decreases gradually (shown in Figure 4(e)). The infill process terminates after the 24th infill point because EI equates to zero for 10 function evaluations successively. It means, it is unlikely to find a new point which improves the current best solution. According to Figures 4(a) and 44(b), the accuracy of the objective and constraint models are both validated, seeing all the sample points lie in the interval [−3,3]. The distribution of sample points is demonstrated in Figure 4(f). Note that only the infill points which have are plotted. As the limited memory LBFGS algorithm will randomly select a point when it cannot find a new point with . As a result, the infill points with are not representative in terms of the objective and constraint values. The *X*-axis is the constraint value while *Y*-axis is the objective value. The dashed red line is the constraint limit. Most of the initial points violate the constraint while only one infill point is beyond the constraint limit, which proves that the infill strategy works well. The point inside the red circle is the optimal solution for the EQ schemes. We solve the EN (Figures 5(a)–55(f)) and RE scheme (Figures 6(a)–66(f)) in the same way.

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

**(f)**

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

**(f)**

**(a)**

**(b)**

**(c)**

**(d)**

**(e)**

**(f)**

##### 3.4. Comparison among Three Schemes

In this section, we compare the results of three schemes. The optimal toll plans for three schemes are illustrated in Figures 7(a)–7(c). The primary *Y*-axis is the distance toll rate for HV groups, while the secondary *Y*-axis is the delay toll rates for AV groups. According to Figure 7(d), all schemes control the network density over time around the critical value (25 veh/km/lane) well. The constraint values under all three schemes are within 2.45 veh/km/lane (shown in Table 3). The RE scheme lets road users to pay 15.79 (10 × 3 AUD) which is 40% higher than the revenues collected in the EQ or EN scheme. However, when we look into the total travel time reduction and the constraint value, the network performance under RE scheme is similar to those under other two schemes. It indicates that RE scheme is less cost-effective. The extremely high toll may induce the low marginal effect.

**(a)**

**(b)**

**(c)**

**(d)**

The EN scheme produces least emissions (7.02 kg). According to the QUARTET pollution emission model, the emission rate goes up with the increase of the cruising speed. Moreover, the optimal toll plan for the EN scheme reduce the emissions by controlling the total time travelling inside the pricing zone, seeing the least total travel time (7726 min.veh) among three schemes. It suggests that the EN scheme can encourage individuals to travel less time in the pricing zone. Nonetheless, since the HVs and AVs are charged differently in our work, in order to increase the acceptance of the policy, should be as small as possible. From this perspective, the EQ scheme outperforms other schemes with AUD/per^{2}. In addition, compared with the no-toll scenario, the EQ scheme reduces the NOx emissions by 18.3% and collects 11.92 (10 × 3 AUD), which is appreciable.

#### 4. Conclusion

The formulation of DCP optimization involves high-dimension decision variables as well as the expensive-to-evaluate objective function. Meanwhile, with the emerging of autonomous technology, HVs and AVs will coexist in the transportation network in the near future. How to design appropriate congestion pricing schemes under mixed traffic is challenging. In this paper, we divide the road users into six groups according to the vehicle type and VOTs. Three congestion pricing schemes with different objectives are investigated. The network performance is considered as a complex constraint to keep the network density close to the critical value. The RK model is adopted to estimate the input-output mapping. Initial and infill sampling strategies are combined to search for optimal toll plans.

A large-scale simulation-based DTA model of Melbourne, Australia, is used to demonstrate the results of proposed pricing schemes. The RE scheme is not advocated because it charges road users too much money without achieving much better network performance. The phenomenon illustrates that demand originates from the pricing zone is evitable. High tolls may obtain low marginal effect. The EN scheme produces least emissions and maintains the relatively high average speed. In this study, AVs are charged by delay while HVs are charged by distance. Meanwhile, people with higher VOTs will have high toll rate. As a result, we should guarantee that the proposed pricing scheme can be accepted by all groups. Hence, the EQ scheme is superior to EN and RE schemes for its smallest . To ensure the efficiency of our pricing schemes.

We should also guarantee that travelers report their VOTs honestly. Nan et al. [15] conducted an in-depth discussion on how to stimulate people to tell truth. For instance, the mechanism design which does not require travelers to provide confidential information. In addition, punishments such as fines or banning the liars from the pricing zone for a certain period are considerable. As long as the government takes a lead role to boost the group-based discrimination, the proposed pricing schemes can be put into practice.

#### Data Availability

The simulation-based traffic model is available from GitHub through https://github.com/meeadsaberi/dynamel.

#### Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this study.

#### Acknowledgments

The first author was supported by the China Scholarship Council grant scheme.