Solving Signal Control Problems with Second-Order Sensitivity Information of Equilibrium Network Flows

The equilibrium network signal control problem is represented as a Stackelberg game. Due to the characteristics of a Stackelberg game, solving the upper-level problem and lower-level problem iteratively cannot be expected to converge to the solution. The reaction function of the lower-level problem is the key information to solve a Stackelberg game. Usually, the reaction function is approximated by the network sensitivity information. This paper firstly presents the general form of the second-order sensitivity formula for equilibrium network flows. The second-order sensitivity information can be applied to the second-order reaction function to solve the network signal control problem efficiently. Finally, this paper also demonstrates two numerical examples that show the computation of second-order sensitivity and the speed of convergence of the nonlinear approximation algorithm.


Introduction
The network signal control problem (NSCP) is to find the optimal signal setting which improves the performance of existing facilities in a transportation network.Conventional methods for optimizing signal settings can be divided into two types: stage-based and group-based approaches [1][2][3][4][5][6].The stage-based approach divided the signal cycle into separate stages and solved the optimal signal settings for each group of compatible traffic movements in stages.This approach is regarded as superior in the concern for safety and loss of capacity with phase switching [6].The group-based approach considered each group of traffic streams having right-of-way in the time domain directly.Compared with the stage-based approach, the group-based approach has a higher degree of flexibility in signal timing arrangement [3,7].However, the most optimization models proposed so far usually converged to a local optimal solution and without taking traffic rerouting effects into account when solving NSCP [5].The equilibrium network signal control problem (ENSCP) is used to find an optimal network signal design when the network flow pattern is constrained to be equilibrium.Friesz [8] points out that this is a problem of interest because of Braess' paradox [9].This paradox shows that the congestion of the network may be severer when adding capacity to a congested network without taking the reaction of network users into consideration.Hence, in practice, the equilibrium network signal design problem must be solved by constraining the network flow pattern to meet user equilibrium.The user equilibrium network design with fixed transportation demand has been studied in both discrete [10] and continuous [11] versions.To help solve the signal control problem, Allsop [12] pointed out that the route choices of road users should be considered the impacts of signal settings changing.Gartner et al. [13] and Fisk [14] described the signal control problem as a Stackelberg or leader-follower game between road users and the administration.The Stackelberg game can be represented as a bilevel problem, where the upper-level problem aims to find the optimal signal setting or link capacity expansions which maximizes system performance, and the lower-level problem aims to solve the user equilibrium (UE) flows, respectively [15,16].
Marcotte [17], Sheffi and Powell [18], Heydecker and Khoo [19], Smith and van Vuren [20], Tan et al. [21], Cantarella et al. [22], Gartner et al. [13], Smith et al. [23], van Vuren and van Vliet [24], and others proposed algorithms to solve the network problem.However, when calculating optimal settings in general road networks, there were no efficient solution algorithms that are combined with anticipating the reactions of road users.Moreover, the iterative optimization assignment (IOA) procedure which solves signal settings and equilibrium flows cannot be expected to converge to the true solution and might lead to a decline in network performance [14,25].The sensitivity analysis-based algorithm evaluates the influence factors as the derivatives of the reaction functions with respect to the upper-level decision variables.The derivative information is obtained by implementing sensitivity analysis for a given solution of the user equilibrium problem [26][27][28][29][30].For the singularity, algorithms can only solve a small network problem.Cho [31] proposed a generalized inverse method, Cho et al. [32] proposed a row reduction method, Patriksson [33], Josefsson, and Patriksson [34] proposed directional derivative method, and Yang and Bell [35] proposed a column reduction method to overcome the singularity problem.In the sensitivity analysis-based linear approximation algorithm, the sensitivity information is used to create a linear approximation of the reaction function and is then inserted into the upper-level problem, iterated until the solutions converge [34,36,37].Recently, Chiou has conducted several studies related to optimal design of area traffic control with equilibrium network flows and proposed a number of computational algorithms to solve the problem, such as the projected Quasi-Newton method and the bundle subgradient projection method [38][39][40][41].Moreover, the ENSCP, combined with an explicit traffic model, TRANSYT, was proposed to evaluate the performance index of the system more precisely [4,38,[42][43][44].
The remainder of this paper is organized as follows: Section 2 introduces the equilibrium network flow models and the first-order sensitivity formula obtained by row reduction method.Section 3 introduces the matrix calculus theory and the second-order sensitivity formula.The network signal control model and solution algorithm are presented in Section 4. Finally, a numerical example and conclusion are presented in Sections 5 and 6, respectively.

Sensitivity Analysis of Equilibrium
Network Flows where ||, ||, and || denote the cardinalities of , , and , respectively.The relationship between arc flow, path flow, and travel demand is given by where Δ is a || × || matrix, with Δ  = 1, if arc  belongs to path  and Δ  = 0 otherwise; Λ is a || × || matrix, with Λ  = 1, if OD pair  belongs to path  and Λ  = 0 otherwise.In general, Δ is called arc/path incidence matrix and Λ is called OD/path incidence matrix.
In sensitivity analysis, a vector of perturbation parameters with dimension ,  ∈   , is introduced.The arc cost function and travel demand function are supposed to be influenced by .Let (, ) be the arc cost function vector and let () be the travel demand function vector.The path cost function vector is given by (ℎ, ) = Δ  (, ).When the network equilibrium is reached, the following equations must be satisfied: where   is the equilibrium path cost of OD pair .Equations (3) are recognized as Wardrop equilibrium conditions; say, there is no traveler can change path unilaterally to improve his travel time [45].Generally, the equilibrium network flow problem can be written in the form of variational inequality (VI) problem as follows [46].Find  ∈ Ω() such that where Ω is the feasible arc flow solution set of the network flow problems.
An equivalent VI can be written with the cost function in terms of path flow variable ℎ rather than arc flow variable  as follows.Find ℎ ∈ Γ() such that (ℎ, )  ( − ℎ) ≥ 0, ∀ ∈ Γ () , where Γ is the feasible path flow solution set of the network flow problems.

First-Order Sensitivity Formula for Equilibrium Network
Flows.The classical first-order sensitivity analysis for equilibrium network flows was proposed by Tobin and Friesz [47].However, Tobin and Friesz method had a strong requirement on the topology of network which may not hold in practical networks [48].Cho et al. proposed the row reduction method to overcome this issue [32].In this paper, we only summarize the key results of the row reduction method, and the readers are encouraged to refer to the original paper [32] for more details.
In the row reduction method, a maximal set of rows from Δ, say Δ 1 , is selected for which the combined matrix [Δ 1 , Λ]  is of full row rank.Hence, Δ can be partitioned as Assume that the number of independent arcs is  1 and the number of dependent arcs is  2 , respectively.Therefore, || =  =  1 +  2 .For a differentiable function,  :  →   , let the partial derivative of  with respect to  (the Jacobian matrix of ) be denoted by   .Let the second-order partial derivative of  with respect to  be denoted by   .The firstorder sensitivity formula can be expressed as where

Second-Order Sensitivity Formula for Equilibrium Network Flows
The sensitivity analysis-based nonlinear approximation heuristic algorithm (NLAA) was firstly proposed by Cho and Lin [37].With the first-order and the second-order sensitivity information, the reaction function of the lower-level problem can be approximated by a nonlinear function.However, Cho and Lin did not provide the general form of the second-order sensitivity information.In this section, we will introduce the preliminary of matrix calculus theory and derive the general form of the second-order sensitivity formula for equilibrium network flows.

Preliminary Definitions and Theorems.
To derive the second-order sensitivity formula, we introduce some definitions and theorems of matrix calculus as follows [49,50].
Definition 1 (Kronecker product).Let  be an  ×  matrix and let  be a  ×  matrix; then the Kronecker product of  and , denoted by  ⊗ , is an  ×  matrix defined by Definition 2 (Vector operator).Let  be an  ×  matrix and   is the th column of ; then vec  is the  × 1 vector: Definition 3. Let  be an  ×  real matrix function of a  ×  matrix of real variables .The derivative of  with respect to  is the  ×  matrix: Theorem 4 (chain rule for matrix functions, Magnus and Neudecker, 1985 [49]).Let  be a subset of  × and assume that  :  →  × is differentiable at an interior point  of .
Let  be a subset of  × such that () ∈  for all  ∈  and assume that  :  →  × is differentiable at an interior point  = () of .Then the composite function  :  →  × defined by () = (()) is differentiable at  and Theorem 5 (Magnus and Neudecker, 1985 [49]).Let  :  →  × and  :  →  × be two matrix functions defined and differentiable on an open set  in  × .Then the simple product  is differentiable on  and the Jacobian matrix is the × matrix: where   and   are the identity matrices of size  and , respectively.
Theorem 6 (Magnus and Neudecker, 1999 [50]).Let  :  →   be a function defined on a set  in   .Let  be an interior point of  and let ( 0 ; ) be an -ball lying in .Let  be a point in   with ‖‖ < , so that  0 +  ∈ ( 0 ; ).If  is twice differentiable at  0 , then the second-order Taylor expansion of function  at  0 +  is where ( 0 ; ) and  2 ( 0 ; ) are the first differential and the second differential of  at  0 , respectively, and 3.2.Second-Order Sensitivity Formula for Equilibrium Network Flows.To derive the second-order sensitivity formula for equilibrium network flows, it is intuitive to take derivative of ( 6) with respect to .For convenience, let where  is an (+ 2 )×(+ 2 ) matrix and  is an (+ 2 )× matrix, respectively.
Lemma 7. The second-order sensitivity for equilibrium network flows is where Proof.Since the first-order sensitivity is the product of ( 17), the second-order sensitivity can be obtained by taking derivative of the product with respect to  directly.According to Theorem 5, the formula of the second-order sensitivity is expressed as ( 16) and the proof is complete.

Network Signal Control Model and Solution Algorithm
In the general problem, the signal variable that can be set by the controlling agent is green time.By specifying the cost functions   for each network arc  in terms of these variables and assuming that the behavioral hypothesis for route choice follows the first principle of Wardrop [45], problem  1 can be presented as If (, ) is strictly monotone, then, for each  ≥ 0, (20) has a unique solution and function () is (continuously) differentiable at every point  ≥ 0. Thus,  2 can be rewritten as  3 : Also, given (, ()) = min ((), ) ⋅ ( − ), then  2 is equivalent to  4 : 4.2.Solution Algorithms.The iterative optimization assignment (IOA) method described by Tan et al. [21] is proceeded as follows.First, fix  and solve (20) for ; then fix  and solve (19) for , continuing this process until  +1 −   → 0 or  +1 −  → 0. The final solution (  ,   ) is termed the Nash solution.Notably, that the solution obtained using the IOA algorithm is not necessarily an optimal solution of the equilibrium network control problem [14].The sensitivity analysis of equilibrium network flows [32,47] was used to solve the equilibrium network signal design problem [26,36,51].

A Sensitivity Analysis-Based Linear Approximation
Heuristic Algorithm.The challenge in solving problem  2 is that, since the lower level of the problem cannot be represented in closed form, it is impossible to obtain an explicit reaction function that can be plugged into the upper level.In the sensitivity analysis-based linear approximation heuristic algorithm, the sensitivity information is used to create a linear approximation of the reaction function and is then inserted into the upper-level problem, iterated until the solutions converge (abbreviated as LAA) [36].
The heuristic is detailed as follows.
Step 0. Determine a fixed small value  > 0 and an initial value  0 .Set  = 0.
Step 4. Reformulate (21) as Step 5. Solve the problem in Step 4 using any software package which can solve the optimal solution for  +1 .If | +1 −   | ≤ , then stop; otherwise set  =  + 1 and go to Step 1.

A Sensitivity Analysis-Based Nonlinear Approximation
Heuristic Algorithm.In the sensitivity analysis-based linear approximation heuristic algorithm, the reaction function of the lower level is based on approximation by a linear function.
In this section, the reaction function of the lower-level problem is based on approximation by a nonlinear function and is plugged into the upper-level problem and is iterated until the solutions converge (abbreviated as NLAA) [37].
Step 0. Determine a fixed small value  > 0 and an initial value  0 .Set  = 0.
Step 3. Using    and   , Taylor expansion and Theorem 6 form the nonlinear approximation  +1 : Since   ,   ,   , and    are known,  +1 can be replaced by a function of  +1 .Thus, Step 4. Reformulate (21) as Step 5. Solve the problem in Step 4 using any software package which can solve the optimal solution for  +1 .If | +1 −   | ≤ ; then stop; otherwise set  =  + 1 and go to Step 1.
In addition to describing the algorithm in more detail, we will provide a proof that if this algorithm converges, it converges to an optimal solution of problem  2 .
So, we know if conditions (i) and (ii) are satisfied, then the following should also be satisfied: (iv)

Numerical Example
This section provides two numerical examples which illustrate the computation of second-order sensitivity and the performance of NLAA.The first example demonstrates the computation of second-order sensitivity in detail.The second example focuses on the speed of convergence between LAA and NLAA.
Example 1.The first example is chosen from Dickson [25] and Fisk [14].The network topology is shown in Figure 1.
The set of OD pairs is {(1, 2), (3,4)} and a signal exists at the intersection of arcs 1 and 3.The cost functions used are where   denotes the green time on arc  and the cycle time,  1 +  3 , is equal to 20.
Additionally, the travel demand  1 from node 1 to node 2 is 10 and the travel demand  2 from node 3 to node 4 is 10.Table 1 lists the arc cost functions   (  ,   ) and the system objective function ().

First-Order and Second-Order Sensitivity Formulas.
In this example,  3 can be replaced by 20- 1 , and  1 that is the only perturbation parameter (control variable) should be considered.Therefore,  is equal to 1. Together with ( 6) and ( 15), the first-order sensitivity with respect to  1 can be rewritten as  where In (36), the sensitivity information of arc flow represents the change of arc flow on arc , respectively, when the control variable  1 increases one unit.Since  1 ≥ 0 and  1 > 0, the equilibrium flow on arc 1 will increase when  1 increases one unit.In the meanwhile, the equilibrium flow on arc 2 will decrease.Because OD pair (3,4) has only one path (arc 3),  1 will not affect the equilibrium flow on arc 3.
From Lemma 7, the second-order sensitivity with respect to control variable is By Theorem 4,    can be derived by the chain rule for matrix functions as follows: In this example, the matrix  is only dependent on  1 .Hence,  ,  = 0 and (38) can be rewritten as where   represents the th column of matrix , and (1 + 2 1 ) , Similarly,    can be derived as According to Definition 1, ( 39), (40), and ( 38), ( 37) can be rewritten as Therefore, 5.2.Computational Results of LAA and NLAA.Based on (36) and (43), the first differential and the second differential of equilibrium arc flow   can be obtained by (14).At iteration , In this example, both LAA and NLAA are implemented in the MATLAB environment.Set  = 0.001 and the initial  1 = 10; Table 2 lists the computational results of LAA and NLAA approaches, and it shows that NLAA is more efficient than LAA. ) and the other is a highway with 5 signal-controlled intersections (arc 1-arc 6-arc 8-arc 10-arc 12-arc 14).The objective of this problem is to find the optimal signal settings which minimize the system cost.
The arc cost functions   (  ,   ) and the system objective function () are listed in Table 3. Table 4 lists the origindestination demand.For the signal-controlled intersections, the arcs entering the same intersection share the same cycle time and the minimum green time for each approach is 10 sec.
In this example, we set  = 0.1 and the initial   =   /2 for each signalized arc.Table 5 lists the computational results of LAA and NLAA, respectively.Two parallel paths from node 1 to node 16 (2-4-16 and 1-6-8-10-12-14) have the same equilibrium travel time 13.4069 min.Table 5 shows that NLAA only takes 6% iterations, compared with LAA, to attain the same level of precision.Figure 3 shows the convergence curves of LAA and NLAA, respectively.The convergence rate of LAA is slower than NLAA due to the zigzag effect.Compared with Example 1, NLAA has more improvement in the speed of convergence than in Example 2. It may imply that NLAA is more efficient to deal with more nonlinear problems.

Conclusions
The key information to solve the equilibrium network signal control problem (ENSCP) is the reaction function of the lower-level problem.Because the reaction function cannot be obtained explicitly, the sensitivity information of equilibrium network flows is used to approximate it.Based on the first-order sensitivity formula and the matrix calculus, this paper first presents the general form of the second-order sensitivity formula for equilibrium network flows.With the second-order sensitivity formula, the reaction function can be approximated more accurately by a nonlinear function.From HSIP to Jhubei city, a simplified real network example demonstrates the speed of convergence between LAA and NLAA.The NLAA has significant improvement in solving the ENSCP with complicated arc cost functions; in this example, the NLAA only takes 6% iterations to attain the same level of precision.This study focuses on the NLAA and a simplified delay formula is adopted to reflect the influence of traffic congestion.Practically, a traffic propagation model, such as TRANSYT model, should be included when solving the ENSCP.Since the derivatives of TRANSYT model have been obtained explicitly [4,42], it can be extended to second-order derivatives and applied to NLAA in the future research.
Compared with LAA, the number of multiplications for matrix multiplication is greatly increasing in NLAA due to the Kronecker-product operation.NLAA has polynomial complexity with the network size and the number of perturbation parameters because of the property of the Kronecker product.There is still opportunity to improve the  computing efficiency through adopting effective Kroneckerproduct algorithms.
respectively.The authors also acknowledge the reviewer's good suggestions for the further research topics.

Figure 1 :
Figure 1: The network topology in Example 1.

Figure 3 :
Figure 3: The convergence curves of LAA and NLAA in Example 2.

Table 2 :
Computational results of LAA and NLAA in Example 1.

Table 5 :
Computational results of LAA and NLAA in Example 2.