Destination and route choice models for bidirectional pedestrian flow based on the social force model

: In this study, social force model (SFM) is extended by using a discretisation grid to permit pedestrians to change their desired speed directions dynamically. In reality, other pedestrians may obscure the visions of the behind pedestrians, so the behind pedestrians will be blocked if they insist to walk in the final desired directions calculated by the SFM. So, a dynamic destination choice model is established to provide pedestrians a series of available intermediate destinations. Based on the dynamic destination choice model, the authors use a discretisation grid to represent all the potential moving directions of pedestrians, and model the weight of every potential moving direction. The direction with the maximum weight is selected as the optimal route at that time step. Besides, pedestrians prefer to pass near to crowds with a low relative velocity and choose the low space occupancy route when several routes have the same pedestrian density. So, pedestrian speeds, space un- occupancy, route length, crowd density and object pedestrian ratio are used to develop the weight model. The modified SFM guarantees pedestrians to obtain an available and optimal route. Compared to other models, the proposed models can be used to reproduce the behavior of bidirectional pedestrians more really.


Introduction
Pedestrian microscopic simulation is an important and increasingly popular method that provides an effective tool for solving some challenging problems. To illustrate the pedestrians' walking mechanism more realistically, various pedestrian microscopic simulation models have been proposed. Social force model (SFM) has been sufficiently validated and widely adopted in pedestrian microscopic simulation.
Most pedestrian microscopic simulation models can be divided into two kinds: (i) models with no route choice and (ii) models with route choice.
In the first type, route choice occurs as self-organisation phenomena which are an emergent phenomenon caused by interaction between pedestrians. These models are not appropriate for complex pedestrian behaviours. On the basis of Newton's second law of dynamics, and considering more realistic behaviours, Helbing and Molnar [1] established initial SFM. SFM quantifies interactions between pedestrians and the outside environment as repulsion and attraction forces, and then updates pedestrians' positions in each simulation time step to simulate the pedestrians' movements. SFM can successfully reproduce many collective pedestrian phenomena such as oscillations at bottlenecks [1], lane formation, clogging at exit doors [2] and faster-is-slower effects [3]. These models, where pedestrians do not have imbedded route choice algorithms, focus on reproducing various observed phenomena in reality. However, when the density of pedestrians in the scenarios is high, simulated results indicated that the simulated pedestrians do not behave as expected, and they move like irrational particles rather than rational pedestrians. Besides, it is common to find some pedestrians heading directly to the final destination set in advance even if the moving direction is occupied by other pedestrians, which results into pedestrians block. In the simulation, pedestrians have no capacity of finding an effective and available route to head to destination. Modelling pedestrians' destination and route choice behaviours are significant tasks in making simulation pedestrian smarter and more realistic.
Pedestrians prefer to choose the shortest and fastest paths, which are the nature of the pedestrians. In most studies, time and distance are the common reasons for choosing a best path. Routes choice can be classified into static and dynamic models. Static route choice models are defined before the simulation begins, and the pedestrians try to walk along the shortest routes. However, pedestrians frequently tend to change their paths according to the complex traffic environment. In dynamic route choice model, the optimal route varies over time. So, dynamic route choice models can represent pedestrian behaviours more effectively. Pedestrians will choose the optimal routes by weighing influencing factors in the environment. Kirik [4] used the floor field model to help pedestrians to find the shortest route by considering the density in front of the moving pedestrians within their sight range, but the strategy of shortest time was neglected. Hoogendoorn [5] adopted the Hamilton-Jacobi-Bellman equation to minimise the travel time. Graph-based routing is a common method used in continuous models (such as SFM), and it is achieved by setting a set of destination points in advance. Geraerts and Overmars [6] proposed a high-quality corridor map method to guide the pedestrian global behaviours and plan routes. The method can obtain a path that is fast and flexible (to avoid conflict with obstacles). On the basis of a graph-based structure, Wagoum et al. [7] proposed an event driven way finding algorithm to represent the pedestrian route choice behaviour. The proposed algorithm helps pedestrians to choose a shortest and quickest route. The shortest route between the start and end positions is unique, but the quickest route changes over time during simulation. So, pedestrians revise their routes according to a travel time-distance analysis. Kretz et al. [8] proposed a route choice model on the basis of minimal remaining travel time. Patil et al. [9] put forward an interactive algorithm to guide pedestrian simulations. The proposed method uses userspecified guidance fields to direct the pedestrians to find a quickest route. Treuille et al. [10] developed a route choice model in which a set of dynamic potentials and velocities were applied to unify collision avoidance and route choice. Asano et al. [11] incorporated a route/exit choice model to give pedestrians the capacity to make far-sighted decisions, making the existing simulation model more realistic. Usually, existing studies tried to obtain the shortest or the fastest routes by optimising the travel time or distance in the constructed network, but few studies tried to develop the models by considering the space occupancy and speed difference. For example, under the same pedestrian density condition, the moving pedestrians may be blocked when other pedestrians stand right in the moving direction. Besides, pedestrians prefer to pass through crowds with a lower-speed difference. How pedestrians deal with these more complex traffic environments is a challenging work. Therefore, more reasonable destination and route choice models are necessary to help pedestrians to find a reasonable and available route rather than only shortest and quickest path. This paper develops a dynamic route choice model by combining distance, density, ratio of object pedestrians and speed difference. We combine modified destination and route choice models into SFM, so that pedestrians can change their desired direction voluntarily at every time step even if the final desired destination is out of sight. Therefore, it is necessary to first establish a dynamic destination choice model to provide pedestrians a series of available intermediate destination. On the basis of the destination choice model, we use a discretisation grid to improve the investigation capacities of the traffic environment, and develop a route choice model based on the discretisation grid. In the established route choice model, pedestrians can evaluate and change their moving directions according to the changing environment. Finally, a comparative analysis between real collected data and simulation data is made. The results show that the established models can be incorporated into the SFM to reproduce more realistic motion.

Social FM
The proposed models are based on the SFM, so we introduce it briefly. The SFM contains four basic forces exerted by other individuals and obstacles: the driving force toward the destination f α 0 (t), the force from obstacles f oα (t), the force from other pedestrians f αβ (t) and random fluctuation force ε α (t). ε α (t) is used to represent diverse random behaviours that cannot be predicted. The resultant force of these forces is applied in a Newtonian equation, driving the pedestrians toward their destinations. The main equations of the SFM are shown in (1) Here, e α (t) is the unit vector of the desired speed direction and e α norm (t) is the vertical unit vector of the desired speed direction. In (2), X ∼ N(0, 1), and its probability density function is as shown in (3)

Driving force
The driving force is defined as a force driving pedestrian α to adopt actual speed v α (t) so as to move in the desired direction e α 0 (t), with a desired speed v α 0 (t) within a relaxation time. It is expressed as (4) where f α 0 (t) is the driving force of pedestrian α at time t, v α 0 (t) is the desired speed of pedestrian α at time t, v α (t) is the current speed vector at time t, τ α is the relaxation time of pedestrian α, e α 0 (t) is the unit vector of the desired speed direction, p α d is the desired destination of pedestrian α and p α c is the current position of pedestrian α.

Interactional force between pedestrians
When a pedestrian approaches other pedestrians, other pedestrians will produce a force against him/her, causing him/her to slow down to avoid collisions with them. The interactional force between pedestrian α and pedestrian β is calculated with (6). The sociopsychological force f αβ soc (t) is used to maintain privacy or personal space from nearby users, and can be expressed as in (7). The physical force f αβ PH (t) contains physical pressure and sliding friction when one pedestrian contacts another, and is expressed as in (9). In general, one pedestrian cannot come into contact with another. If physical force occurs, an accident will happen where f αβ (t) is the force between pedestrians α and β, A αβ and B αβ are the strength coefficients to be estimated, p α (t) and p β (t) are the centre coordinates of pedestrians α and β, r αβ is the sum of the radiuses of pedestrians α and β, d αβ is the distance between the centres of pedestrians α and β, n αβ is the vector pointing from pedestrian α to pedestrian β, F αβ is the anisotropic factor, Δv αβ t is the speed difference between pedestrians α and β in the tangent direction, t αβ is the vector in the tangent direction and K and k are very large constants.
F αβ was introduced by Helbing and Molnar [1], and it can be calculated as shown in (13). When λ α takes a larger value in the range 0 ≤ λ α < 1, it means that the pedestrian has a larger influence on the pedestrian in front of him/her than the pedestrian behind

Interactional force between pedestrian and obstacle
The interactional force between pedestrians and obstacles also contains the socio-psychological force f oα soc (t) and the physical force f oα PH (t). The interactional force between pedestrian α and an obstacle is calculated as shown in (16)  Here, n oα is the vector to the centre of an obstacle and t oα is the vector in the tangent direction.

Route choice model
Pedestrian' position at every simulation time step is determined by the magnitude and direction of speed. The destination model provides a series of optional destinations, and one optional destination represents one desired speed direction. It means that pedestrians face many optional desired speed directions. How to choose an optimum desired speed direction is a route choice problem. As shown in Fig. 2a, a discretisation grid helps every pedestrian to find the optimal direction of movement and most efficient route. Inspired by Guo and Huang [12], the space around the pedestrian is divided into n possible directions with different weights. In every simulation time step, the direction with the largest weight is selected as the desired speed direction.
In this paper, pedestrians can re-evaluate and change their intended directions in each iteration. To choose the most efficient route to the next intermediate desired destination, the space around each pedestrian is divided into a discretisation grid. On the basis of the studies of Guo and Huang [12] and Saboia and Goldenstein [13], we model the weights of each grid, and propose a new model of occupation of the space surrounding each pedestrian.
As shown in Fig. 2b, each pedestrian i occupies a real space that can be represented by a circle with radius R i . Here, η is the minimum distance that a pedestrian expects to keep between him/her and other individuals. Besides this, each pedestrian i has a virtual space, which represents the pedestrian's greatest neighbourhood, denoted by PR. The annulus area is formed by an inner circle with radius R i + η and an external circle with radius PR. The discretisation grid is the annulus area, and the annulus area is divided into n lattices. c k is used to denote each lattice, where k ∈ 0, n − 1 . Fig. 3a represents a discretisation grid of a pedestrian space. L 0 , L 1 , …, L n − 1 is a series of counterclockwisedirection sequence boundary vectors. One lattice has two boundaries, L a and L b . The opening angle of lattice c k is θ = 360 ∘ /n. The direction of lattice c k is denoted by the vector d i k .
Taking the vector d i s [calculated in the same way as shown in (5) (Fig. 3c).

Weights of discretisation grid
For a pedestrian i, the weight of each lattice c k is denoted by w i k .
Here, w i k can be calculated with (20) The term N is a normalisation coefficient that is responsible for guaranteeing that ∑ k w i k = 1. The term δ i k represents people's preference for directions in which there are fewer object pedestrians, lower densities and a shorter route to the intermediate destination. It can be estimated using (21) where d route represents the length of the shortest route from the current position of pedestrian i to one intermediate destination, ρ is the pedestrian density in lattice c k , Area Lattice c k represents the area of lattice c k , p object represents the proportion of object pedestrians out of all pedestrians in lattice c k , n object represents the number of object pedestrians in lattice c k and n subject represents the number of subject pedestrians in lattice c k . Considering that pedestrians prefer to pass near to crowds with a low relative velocity, E i k is introduced to represent this characteristic. E i k is evaluated as the difference between the current velocity v α (t) of pedestrian α and that of other pedestrians β, v β (t), walking in the lattice c k . Only pedestrians within the lattice c k are taken into account in estimating the speed difference. Pedestrian α prefers to pass through crowds with a lower-speed difference. E i k can be estimated using (25) In (20), the term F i k is similar to the anisotropic factor. It is used to represent the strength that makes one pedestrian take the direction d i k , and is a drift coefficient varying between 0 and 1.
The term Access i k represents the strength of the lattice c k that is being made use of by pedestrian i, for example, how many pedestrians are occupying the path of pedestrian i in lattice c k . Access i k is expressed as in (27) Access where w i is all the pedestrians in lattice c k , and s i j k is the strength of pedestrian i departing from j in lattice c k , takes values from 0 to 1, and can be expressed as in (28) where A i, j k is the ratio of separation between pedestrian i and the nearest pedestrian j. Here, j must have an intersection with lattice c k . A i, j k can be calculated using (29) where H i, j k is the distance between pedestrian i and the nearest pedestrian j in lattice c k , as shown in Fig. 4b. In (28), the term B i, j k is the space occupancy of lattice c k . It is the angular opening of the lattice c k , where the path of pedestrian i is not occupied by pedestrians lying in the lattice c k .
Here, we take one lattice c k as an example to illustrate the process of calculating B i, j k of lattice c k in each iteration. As shown in Fig. 4c, a vector sweeps the area of lattice c k , from boundary L k to boundary L k + 1 in a counterclockwise direction. When the vector meets the first pedestrian in lattice c k , the beginning tangent line of the first pedestrian is denoted by l 1 ; then, the vector goes on sweeping the remaining area of lattice c k . If there are no pedestrians in the direction of the end tangent line of the first pedestrian, the end tangent line of the first pedestrian is denoted by l 2 . Otherwise, the vector goes on sweeping the area of lattice c k until a tangent line of pedestrian is not occupied by pedestrians. For example, as is shown in Fig. 3c, the sweeping vector meets the second pedestrian, the beginning tangent line of the second pedestrian is denoted by l 3 , and then the vector goes on sweeping the remaining area of lattice c k , because the end tangent line of the second pedestrian is occupied by the third pedestrian and the fourth pedestrian, so the end tangent line of the second pedestrian is not recorded. The vector goes on sweeping the remaining area of lattice c k , because the end tangent line of the third pedestrian is occupied by the fourth pedestrian, so the end tangent line of the second pedestrian is also not recorded. The vector goes on sweeping the remaining area of lattice c k , and the end tangent line of the third pedestrian is not occupied by pedestrians, so the end tangent line of the fourth pedestrian is denoted by l 4 . The vector goes on sweeping the remaining area of lattice c k until it meets the boundary L k + 1 . After sweeping the area of lattice c k , a series of vectors l 1 , l 2 , …, l n is obtained. The space between vector l 1 and vector L k (vector l 3 and vector l 2 ; vector l 5 and vector l 4 ; vector l 6 and vector L k + 1 ) can be made use of by pedestrian i to head for one intermediate destination. The opening angle α 1, k (α 2, 3 , α 4, 5 , α 6, k + 1 ) between vector l 1 and vector L k (vector l 3 and vector l 2 ; vector l 5 and vector l 4 ; vector l 6 and vector L k + 1 ) is calculated. The bisector of the two vectors with maximum opening angle of α 1, k , α 2, 3 , α 4, 5 , α 6

Decision making over new desired direction
The direction d i k of lattice c k with the maximum weight will likely be selected by pedestrian i as his/her new potential path e new, i for the next simulation time step. The desired velocity is expressed as in (32). When the weight of the new route is larger than that of the current route of a pedestrian with worthwhile constant g i , the pedestrian will make a decision to change his/her current route to the new route. The decision-making model is expressed as in (33) where e current, i 0 (t) is the unit vector of the current speed direction; w current, i k is the weight of the current moving direction; w new, i 0 is the weight of the new moving direction.

Conflict avoidance
Asano et al. [11] assumed that the subject pedestrian would adopt a giving-way manoeuvre to avoid collisions with the nearest pedestrian if they stepped into his/her private sphere. Two common conflicts are shown in Fig. 5a. However, since not all conflicts turn into collisions, we need to judge whether the conflicts are valid. As described in the study by Zeng et al. [14], the time to conflict point (TTCP) can be used to identify whether a conflict is valid or not. The TTCP represents the time it will take for two pedestrians to pass the CP at their current speeds and directions.
The TTCPs for the subject pedestrian α and the conflicting pedestrian β are calculated using (34) and (35). Only when the TTCPs of both pedestrians are positive, the conflict between them can be identified as valid. A negative TTCP value signifies that a pedestrian has passed the CP and no collision will occur  T αβ is the absolute value of the time difference between now and the CP, and is expressed as shown in (37); T αβ can be used to evaluate the collision risk between two pedestrians. Larger values of T αβ signify lower collision risk Subject pedestrian α will adopt a giving-way manoeuvre to avoid valid collision with conflicting pedestrian β. As shown in Fig. 5b, the subject pedestrian takes measures to avoid collisions by adjusting his/her speed or direction. The giving-way manoeuvre guarantees that two pedestrians can go on moving after giving way, rather than stopping forever due to force balance, which needs to be avoided. In the real world, the giving-way manoeuvres adopted by pedestrians are flexible and random. The pedestrians can pass by the CP by slowing down or speeding up. At the same time, whether or not a pedestrian changes direction is dependent on their random behaviours. Zeng [14] proposed a stochastic adjustment mechanism to reflect direction changing. This paper also adopts this manoeuvre.
The adjusted speed vector v α * can be expressed by (38) where ψ is the anticlockwise rotation angle and Δv is the speed increment. Therefore, the giving-way manoeuvre can be presented as an individual force f αβ r from the nearest conflicting pedestrian, and can be calculated as in (40)

Simulation results and discussion
Since there are many parameters contained in the proposed models, and these parameters determine the accuracy of simulation model. So, model parameters should be calibrated with observed data first.

Data acquisition
The East crosswalk at Tongzhi Street in Changchun City, China was selected as the study site. By applying the image-processing software developed by Jiang [15], the crossing coordinate positions of bidirectional pedestrians were manually obtained from the video every 0.5 s. In total, the crossing trajectories of 542 pedestrians were obtained and 57,600 crossing coordinate positions of those bidirectional pedestrians. Software interface for data acquisition is shown in Fig. 6. Velocities and acceleration can be obtained directly. Distances and velocity angles between pedestrians can be derived by calculation.

Calibration of model parameters
There are many parameters in the modified SFM. Some of the parameters can be obtained from the observed data directly, but others cannot. The maximum log-likelihood estimation method was used by Zeng et al. [14] to calibrate the parameters of modified SFMs. The calibration method used in this paper is based on the method adopted by Zeng et al.
In the study by Zeng et al. [14], a three-stage process is adopted to overcome this problem: (i) The parameters K, k, λ α , r i , PR and m, which are measurable, but difficult to identify from the trajectory dataset, are estimated by reference to the research results of Zeng et al. [14], Anvari et al. [16].
(ii) Those parameters that can be determined from the observed dataset, for example, relaxation time τ α and desired speed v α 0 (t), are estimated directly.
v α 0 (t) is the desired speed, and it takes the value of pedestrians' average speed when they cross the crosswalk with no other pedestrians or turning vehicles during the pedestrian green light. As shown in Fig. 7, the average desired speed is 1.37 m/s and the standard deviation is 0.27. (iii) After determining these measurable parameters, nonmeasurable parameters such as the strength coefficients of each force that do not have concrete physical meanings are derived by the maximum log-likelihood estimation method.
Real trajectories are extracted from video data. The trajectory extractor software developed by Jiang was used to extract the trajectories, velocities and accelerations of pedestrians from video recordings. The error and accuracy of the software were analysed in detail, and the video analysis procedure was explained by Jiang [15].
Other non-measurable parameters, namely A αβ , B αβ , A oα , B oα , F, η and g i , are calibrated by the maximum log-likelihood estimation method. The next time step position [P α (t k + 1 | θ p )] is determined by the resultant force, and the resultant force closely depends on the SFM parameters (θ p , which represents the set of all model parameters, A αβ , B αβ , A oα and B oα ). The distance vector Δd α est (t k | θ p ) points from P α (t k | θ p ) to P α (t k + 1 | θ p ). Here, the x and y components of Δd α est (t k | θ p ) are assumed to obey a normal distribution with mean (μ x , μ y ) and standard deviation (σ x , σ y ). (μ x , μ y ) and (σ x , σ y ) can be calculated from the observed Δd α obs (t k | θ p ). In a single prediction step (from t k to t k + 1 ), the likelihood L k is expressed as in (41) Since the x and y direction components are assumed to obey normal distributions, Δd α (θ p ) can be expressed as in (42) μ and σ can be calculated using (43) and (44) where Δd αx est (t k | θ p ) is the x component of Δd α est (t k | θ p ) and Δd αy est (t k | θ p ) is the y component of Δd α est (t k | θ p ).
The L k (Δd α (θ p ) | μ, σ) of model parameters θ p are obtained such that (44) is maximised The calibration results of the parameter estimation are presented in Table 1. These parameters are applied in simulating the developed SFM. The p values at the 95% confidence level indicate that all the parameters of the model are effective.

Example analysis
A simulation analysis of bidirectional pedestrian flow was carried out. The length of the corridor is 40 m. The width of the corridor is 8 m. To acquire stable simulation results, the simulation was repeated 100 times to eliminate the influence of random disturbances. Pedestrians enter the corridor from both sides at a rate of ten pedestrians per second. All the pedestrians access the corridor with the same desired velocity of 1.37 m/s. The anisotropic character λ α of the pedestrians is set to 0.3, so that interactions outside of the field of view have little effect on them.
After 100 simulations were performed, all the average speeds of individuals in the simulation were recorded. The speeds of all pedestrians in the 100 simulations were used to determine the speed distribution. Then, the simulation speed distribution and the observed speed distribution were compared. Fig. 8a shows the analysis of the pedestrians' speed distribution. It is clear from the comparison of their means and standard deviations that the simulation and observed speed distributions are closely matched.
As shown in Fig. 8b, the observed data (velocity-density points) were compared with the simulation data from Helbing's model [1] and our model. It is obvious that the simulation data from our model is in greater agreement with the collected data. The simulation speeds from Helbing's model [1] are smaller than the collected data. However, with a slight reduction in speed, particularly when the pedestrian density is <1 (pedestrian/m 2 ), the difference between the simulation data from Helbing's model and our model is very small. Since pedestrians have good vision and enough space to walk in low densities, the route choice model has little meaning. However, as density increases, especially when it is larger than 1 (pedestrian/m 2 ), the route choice model will improve the pedestrians' walking speed. The higher the density, the more benefits are gained from the route choice model. The results are consistent with the realities of bidirectional pedestrian flow. Fig. 8c shows the evacuation time from the first pedestrian stepping into the simulation corridor to the last pedestrian stepping out of the simulation corridor. It is clear that, when the number of pedestrians is constant, the evacuation time of pedestrians in Helbing's model is larger than that in our model. However, when the number of pedestrians is small, the evacuation times of the two models show little difference, which is due to the pedestrians having good vision and enough space to walk in low densities when the number of pedestrians is small and the interference between them is weak. This tendency is similar to the tendency of the speed-density diagram. Only when the number of pedestrians increases past a certain amount, do the pedestrians start to interfere with each other. As the number of pedestrians increases, the evacuation times of the bidirectional pedestrians will increase.
Since the desired destination in Helbing's model was set as a single point which the pedestrian would head for. It results in many unreasonable behaviours in the simulation. Pedestrian will stop when his path is blocked rather than choosing another path to move, which leads to lower speed and more evacuation time. Pedestrians in our models, pedestrian can change his desired destination and route at every time step, so they can move at a higher speed most of the time with less stop.
We can see the potential application of the developed model from Fig. 8c. For example, we can draw the conclusion that the route choice model is effective only when the number of pedestrians exceeds 54 or the local density is larger than 1. This is because the pedestrians scarcely interact and have little influence  on the evacuation times of other pedestrians when there are <54 of them in the simulation corridor.

Conclusions
Simple local models are significant for improving the decisionmaking capacity of individual pedestrians in a simulation. In this paper, two improved tactical aspects of destination and route choice models have been combined into the SFM. This paper focuses on the destinations, route choice and conflict avoidance of pedestrians. The pedestrians use their vision field to avoid counter flow and consequently to find the most efficient route to take. Before building the route choice model, we developed a dynamic destination choice model to assign pedestrians a series of intermediate destinations. Then, we adopted a discretisation grid to permit the pedestrians to change their intended direction dynamically from time step to time step based on the evaluation of the environment around them. Several validations of the modified SFM were conducted by comparing the simulated and observed speed distributions and analysing the relationship between speed and density. The comparison shows that the modified SFM has three merits: first, the bidirectional pedestrian flow with the route choice model is more efficient; second, anticipated blockages are less likely to happen; and third, our model can deal with conflict avoidance effectively. These results are consistent with the real data on bidirectional pedestrian flow. Future study will focus on extending our models to simulate pedestrians in complex environments such as markets and subway stations.