Implementation PT Demand - smart-fm/simmobility-prod GitHub Wiki
1 References
[1]“Tech_Specs_Public_Transport_Demand_DynaMIT_SimMobility”, version 0.4
[2]Yen, J. Y. (1971). Finding the k shortest loopless paths in a network. management Science, 17(11), 712-716.
[3]Tan, Rui; Robinson, S; Lee, D.H.; Ben-Akiva, M (2014) Evaluation of Choice Set Generation Algorithms for Modeling Route Choice with Smart Card Data in Large-Scale Public Transport Network.
2 Network Data Model
2.1 Input Data
The Network data used is from Google Transit Network Data. It consists information on stops and service line schedules including trains (MRT, LRT) and buses.
Train stops are having the following format as listed in the following table.
| Variable | Description |
|---|---|
| Stop_id | stop id for stops as specified by LTA. Examples: “NS27”, “CC15/NS17” |
| Stop_code | stop code for stops as specified by LTA. Examples: “NS27”, “CC15/NS17” |
| Stop_name | Name of stops as specified by LTA. Examples: “Marina Bay” |
| Stop_desc | Description of stops, usually the street where the stop is located. Example: “Marina Street” |
| Stop_lat | Latitude of stops. Example: 1.2761 |
| Stop_lon | Longitude of stops. Example: 103.85467 |
| Zone_id | MTZ zone where the stop is in, as specified by URA.,Example “MRB” |
Bus stops are having the following format as listed in the table below.
| Variable | Description |
|---|---|
| Stop_id | stop id for stops as specified by LTA. Examples: “14429”, |
| Stop_code | stop code for stops as specified by LTA. Examples: “14429” |
| Stop_name | Name of stops as specified by LTA. Examples: “CP B” |
| Stop_lat | Latitude of stops. Example: 1.263798333 |
| Stop_lon | Longitude of stops. Example: 103.8041983 |
Train service schedule information and bus services schedule information is extracted from the following data:
| Variable | Description |
|---|---|
| trip_id | Service trip id. Example “100_saturday_1-S” |
| Arrival_time | Scheduled arrival time at the specified stop in stop_id. Example: “05:30:00” |
| Departure_time | Scheduled departure time at the specified stop in stop_id. Example: “05:30:00” |
| Stop_id | Stop id for the current stop along service line as specified in trip_id. Example “66009” |
| Stop_sequence | The sequence number for current stop as specified in stop_id along the service line as specified in trip_id |
2.2 Build Network
The final vertex and edges for network built up are of the following format respectively:
| Variable | Description |
|---|---|
| Stop_id | stop id for stops as specified by LTA. Examples: “14429”, “NS27” |
| Stop_code | stop code for stops as specified by LTA. Examples: “14429”, “NS27” |
| Stop_name | Name of stops as specified by LTA. Examples: “CP B”, “Marian Bay” |
| Stop_lat | Latitude of stops. Exmaple: 1.263798333 |
| Stop_lon | Longitude of stops. Example: 103.8041983 |
| EZLink_Name | The stop name as recorded in EZlink Card data. Example, “14429” “STN Marina Bay” |
The final edges are having the following format:
| Variable | Description |
|---|---|
| Start_stops | Stop id for the starting vertex |
| End_stops | Stop_id for the ending vertex |
| R_type | Service line type. Example: “Bus”, “LRT”, “Walk” |
| Road_Index | Index for road type: 0 for Bus, 1 for LRT, 2 for Walk |
| Road_edge_id | Strings of passing road segments by current edge as specified in edge_id. Example:”4/15/35/43”, each number corresponds to a road segment between two adjacent stops along service lines. Road_edge_id containing single “edge_id” are referring to the road segment between two adjacent stops along service lines as specified in R_service_lines. Road_edge_id containing more than 1 “edge_id” are virtual edges representing traversing the edges sequentially without a transfer in between |
| R_service_lines | if the edge is a route segement, it contains services line/direction that travelling on this route segment. If this edge is a walking leg, it is string “Walk” |
| Edge_id | Id for current edge |
| Link_travel_time | Estimated travel time on current edge. If it is a transit leg such as RTS or Bus, it is the estimated in-vehicle travel time on current route segment. If it is a walking leg, it is the estimated walking time on current walking leg. |
| Transfer_penalty | Transfer penalty-used to impose penalty on transfers for shortest path searching which assumes path cost is the addition of edge attributes. |
| Waiting_time | Estimated waiting time to embark on current edge. Data obtained from scheduled service line information. |
Note that the edge here is corresponding to both the route segments and walking legs.
iGraph package in R is utilized to build up the graph for shortest path searching with input of vertex and edges as specified in tables. An embedded shortest path searching algorithm is used for shortest path searching in which Bellman-Ford algorithm is used to allow parallel edges between two vertices. The output path format can be either a list consisting of edge ids in sequence or vertices in sequence. Due to the existing of parallel edges in the network, the path format is restricted to list of edge ids in this work.
Below is an example of R code of building network based on vertices and edges provided, as well as searching for shortest path.
- #build up network based on vertices (stops) and edges (links)
- temp_g<-graph.data.frame(links, directed = TRUE, vertices = stops)
- #searching shortest path from origin to dest with edge attributes specified in weights. Output = “epath” specifies the path format is a list of edges along the path.
- new_path<- get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty+E(temp_g)$waiting_time, output= "epath")
3 Implement Stop Pair Choice Set
This section presents how to generate stop pair choice set given an origin and destination. Here, the stop pair means an initial boarding stop and a final alighting stop. The origin and destination is given from pre-day module with exact geometry location in terms of latitude and longitude. Firstly, aA selection criteria based on geometry regions will be adopted for stop pair selection.
3.1 Geometry Selection of Stop-Pair
Denote the geometry points of given OD pair is characterized by
- origin$lat for latitude of origin point
- origin$lon for longitude of origin point
- dest$lat for latitude of destination point
- dest$lon for longitude of destination point
3.1.1 Selection Criteria
Denote the geometry location of stops as stop$lat, and stop$lon, Any SStops (Bus stop or MRT/LRT stations) that are within satisfies the following condition is selected as a candidate boarding stop with respect to the origin:
- Origin$lat-g_interval<=Stop$lat<=origin$lat + g_interval
- Origin$lon – g_interval <= Sstop$lon <= origin$lon + g_interval
Any Stop (Bus stop or MRT/LRT stations) that satisfies the following condition is selected as a candidate alighting stop with respect to the destination: - Dest$lat-g_interval<=Stop$lat<=dest$lat + g_interval
- Dest$lon – g_interval <= Sstop$lon <= dest$lon + g_interval
All combinations of candidate boarding stop and candidate alighting stop forms the stop pair choice set.
3.1.2 Pseudo code
Given OD, origin$lat, dest$lat, origin$lon, dest$lon,
Let g_interval = 0.02 ( can be adjusted accordingly, 0.02 is to select a square with 500-meter radius)
Initialization:
g_interval = 0.02
Boarding Stop Set BS={}
Alighting Stop Set AS={}
Stop pair SP={}
Search Boarding Stop Set BS:
Let Search through all vertexes,
If stop I satisfy:
- Origin$lat-g_interval<=Stop$lat<=origin$lat + g_interval
- Origin$lon – g_interval <= Stop$lon <= origin$lon + g_interval
Then
BS<-BS union {stop}:
Search Alighting Stop Set AS:
Search through all vertexes,
If stop I satisfy:
- Dest$lat-g_interval<=Stop$lat<=dest$lat + g_interval
- Dest$lon – g_interval <= Stop$lon <= dest$lon + g_interval
Then
AS <-AS union {stop}
From Stop-Pair Choice Set for OD:
For each Boarding stop B_i in BS
For each alighting stop A_j in AS
If B_i!=A_j,
SP<-SP union {B_i, A_j}
Compute access walk from origin to B_i and egress walk from A_j to destination
Store them as attribute with this stop pair.
BS<-BS union {stop}
4 Implementing stop-to-stop Choice set generation methods
Four choice set generation methods will be implemented, and there are: Labeling approach, Link Elimination Approach, K shortest Paths Approach and Simulation Approach.
Path Feasibility
A path feasibility check is imposed on all approaches to filter out unfeasible paths. The feasibility check ensures that for each path, the first link and the last link should not be a walking link and only a single walking link or a direct transfer is allowed between two stops.
4.1 Desired Path format
As the shortest paths generated at this step will be firstly concatenated with an origin point and destination point, and subsequently send to supply side for simulation. The desired path format is a list of transfer points along the path.
The desired stop-to-stop paths choice set PSsSt for stop-to-stop pair {Ss,St} contains set of paths specified in the following format:
Where Pi is the ith path in PSsSt, Ss is the starting stop, S1i,S2i,S3i,....,Sqki are the transfer nodes in sequence where Qk denotes the number of the last transfer node before reaching ending stop St. A necessary conversion from the generated path format by embedded shortest path search algorithm to the desired format is needed.
4.2 Labeling Approach
Labeling approach searches shortest paths based on different cost functions as defined as “Labels”. In current implementation, the following labels can be included.
| Labels | Label Description |
|---|---|
| 1 | Minimal total in-vehicle travel time |
| 2 | Minimal number of transfers |
| 3 | Minimal walking distance |
| 4 | Maximized travel on MRT |
| 5 | Maximized travel on Bus |
| 6 | Minimal waiting time at transfers |
| 7 | Minimal total travel time1 (in-vehicle +waiting +walking) |
| 8 | Minimal total travel time 2( in-vehilce+5waiting+3walking) |
| 9 | Minimal total travel time 3( in-vehilce+5waiting+5walking) |
| 10 | Minimal total travel time 3( in-vehilce+5waiting+8walking) |
4.3.1 Pseudo code
Below is an example of pseudo code for labeling approach:
Initialization:
Searched Paths Q=empty.
Search Index i=1
Procedure:
Step 1: Generate cost function of label i,
Step 2: Search shortest path p_i based on current cost function
Step 3: If p_i does not exist in Q, Q=Union{Q, p_i)
Step 4: Increment I and check whether stop criteria has been met
Stop Criteria:
I>10
Cost function for each Label can be configured as follows with respect to the link attributes.
Label 1. Minimal total in-vehicle travel time:
- a. Set Link_travel_time =0 for Walking Legs (R_type==”Walk”
- b. Set transfer penalty= C (constant, non-zero) for all
- c. Cost function = (link_travel_time+transfer_penalty)
Label 2. Minimal number of transfers - a. Set Transfer_panalty = C (constant, non-zero) for all
- b. Cost function = transfer_penalty
Label 3. Minimal walking distance - a. Set Link_travel_time= C (constant, non-zero) for non-Walking Legs(R_type!=”Walk”)
- b. Set transfer penalty= zero for all
- c. Cost function = Link_travel_time
Label 4. Maximized travel on MRT - a. Set Link_Travel_time = B (constant, big value) for all bus route segments (R_type==”Bus”)
- b. Cost function = Link_Travel_time+transfer_penalty+Waiting_time
Label 5. Maximized travel on Bus - a. Set Link_Travel_time = B (constant, big value) for all RTS route segments (R_type==”RTS”)
- b. Cost function = Link_Travel_time+transfer_penalty+Waiting_time
Label 6. Minimal waiting time at transfers - a. Cost function = Waiting_time
Label 7. Minimal total travel time1 (in-vehicle +waiting +walking) - a. Cost function = Link_Travel_time+transfer_penalty+Waiting_time
Label 8. Minimal total travel time2 (in-vehicle +5waiting +3walking) - a. Set Link_Travel_time=3*link_travel_time for walking legs (R_type==”Walk”)
- b. Set waiting_time=5*waiting_time
- c. Cost function = Link_Travel_time+transfer_penalty+Waiting_time
Label 9. Minimal total travel time3 (in-vehicle +5waiting +5walking) - a. Set Link_Travel_time=5*link_travel_time for walking legs (R_type==”Walk”)
- b. Set waiting_time=5*waiting_time
- c. Cost function = Link_Travel_time+transfer_penalty+Waiting
Label 10. Minimal total travel time4 (in-vehicle +5waiting +8walking) - a. Set Link_Travel_time=8*link_travel_time for walking legs (R_type==”Walk”)
- b. Set waiting_time=5*waiting_time
- c. Cost function = Link_Travel_time+transfer_penalty+Waiting
4.4 Link Elimination Approach
The main purpose of link elimination is to demonstrate how a passenger will react when there is a certain blockage in the network. Link elimination is analogous to simulating a blockage at certain road segments. During the search for the least cost path it is assumed that people will always chose the current least cost path given the current network conditions.
Note that for Link Elimination Approach, the filtering criteria cannot be simply applied when all the paths have been found. Making sure the first shortest path is feasible is very important for Link Elimination method as the searching of subsequent paths is dependent on the first shortest path. Therefore, it is recommended that some actions should be taken to ensure that the first shortest path found is feasible, meaning it satisfies the feasibility criteria.
The following two data tables are needed for Link Elimination Approach: 1) Match_Table_EtoR and 2) Match_Table_RtoE.
| Variable | Description |
|---|---|
| edge_id | Id for current edge, the same as in edge data table |
| road_edge_id | Strings of passing road segments by current edge as specified in edge_id. Example:”4/15/35/43”, each number corresponds to a road segment between two adjacent stops along service lines. Road_edge_id containing single “edge_id” are referring to the road segment between two adjacent stops along service lines as specified in R_service_lines. Road_edge_id containing more than 1 “edge_id” are virtual edges representing traversing the edges sequentially without a transfer in between. |
Match_Table_RtoE identifies the actual road segments by the specified edge_id, while Match_Table_EtoR tabulates the affected edge ids in strings given there is a blockage happens at the road segment as specified in road_id.
| Variable | Description |
|---|---|
| road_id | Id for road segments, as listed in road_edge_id in edge data table and Match_Table_EtoR table.,Each Id refers to an actual road segment between two adjacent stops, which is traversed by at least one service line/direction. |
| affected_edge_id | Strings of affected edge ids given a blockage happens at the road segment as specified in road_id.,Example:” 1/6670/6674/6675/6676/6677/6678/6679/6680”, each number corresponds to an edge in edge data table. |
4.4.1 Pseudo code
Below is an example of pseudo code for Link Elimination Approach:
Initialization:
Construct Network g_start from vertex data and edge data.
Searched Paths Q=empty.
Get path origin stop : Origin, destination stop: Dest
Search Index i=1
Eliminate index j=1
Procedure:
Step 1. Search for a feasible least cost path denoted as path p_1
One possible way of ensuring feasibility:
a. Set Link_Travel_time = B (constant, big value) for all walking legs with {R_type==”Walk” & (origin == Origin | dest==Dest)}
b. Search for shortest path p_1
c. If p_1 ==NULL, then STOP, otherwise, Q={p_1} , get total number of edges along p_1: J, denote the network as g_start, go to Step 2
Step 2. Resume network to g_start,
Step 3. Identify the associated road segment along least cost path p_1, and its affected edge ids if a blockage happens on the road segment E.
a. Get the j-th edge id along p_1, denote e_j
b. if e_j is a walking link, denote affected_edge_ids E={e_j}
c. If e_j is not a walking link, get the road_edge_ids from Match_Table_EtoR, denote R, for each r in R, get the affected_edge_ids from Matched_Table_RtoE, denote as E
Step 4. Simulate road blockage on associated road segment by setting the link_travel_time = B (constant, big value) to all edges in E in the network.
Step 5. Research for least cost path p_i under current network.
Step 6. If p_i is feasible and does not exist in Q, Q=Union{Q, p_i), increment i, increment j.
Step 7. If stop criterion is met, then stop. Otherwise, go back to Step 2.
Stop Criteria:
All link along p_1 has been eliminated once ( j>J)
4.5 K shortest Path Approach
The K-shortest path algorithm of reference 2- Yen (1971) is adopted in this paper, with the following modifications to make it suitable for multimodal PT network: 1) instead of comparing node sequences, it has been modified to compare link sequences for determining links to eliminate, as parallel links are existing in our PT network; 2)all walking links at origin are eliminated. K is set to 30 in this paper as the maximum size of observed choice sets is only 15 from the smart card data.
4.5.1 Pseudo code
Initialization:
Cost function = Link_Travel_time+transfer_penalty+Waiting_time
Search shortest path p_1,
Searched Path Q={p_1}
temp paths list B={empty}
Number of Paths limit K,
Search Index k=1
Procedure:
Step 1:Increment k-Search for k-th shortest path:
Identify edge sequences to compare from p_(k-1). Denote edges in p_(k-1) as edge sequence{e_1, e_2, e_3,….e_q(k-1)} where e_q(k-1) is the SECOND last edge in p_(k-1) in sequence. q(k-1) is the [number of edges in p_(k-1) -1]
- Initial search
a) Eliminate first edge as in {p_1, p_2,…p_(k-1)} from the network
b) Search for shortest path, denote as p_temp, if p_temp does not exist in Q, or B, B=union{B, p_temp}. Otherwise, do nothing.
- If q(k-1)>=1, start subsequent search:
For i in 1: q(k-1)
{
a) check whether the sub-path {e_1, e_2,…,e_i} coincides with sub-path consisting of the first i edges in sequence for j =1 ,2,3, k-1.If coincides, eliminate edge e_(i+1) along that path.
b) Denote subpath {e_1, e_2, …e_(i-1)} as S_i, Search for shortest path R_i with origin at Ending Vertex of (edge_i) and destine at the same destination.
c) concatenate S_i and R_i to obtain A_i_k,
if A_i_k does not exist in B nor exist in Q,
add A_i_k to B
B=union{B, A_i_k}
}
}
Step 2: Compute path cost for each path in B, record the least cost path A_j as p_k.
Remove A_j from B
Stop Criteria:
B is empty, or k=K
4.6 Simulation Approach
In-vehicle travel time and walking time are both randomized following an independent and identically distributed normal distribution with mean equals to its original value and standard deviation is set to 5 times the original value. To avoid drawing a negative value, the absolute value is taken. 50 draws of randomized travel times were performed for each sample OD. The selection of sampling distribution and number of draws takes into consideration of the maximum size of observed choice set, coverage, and computational time.
4.6.1 Pseudo code
Initialization:
Draws = N, assign distribution = f(u, m)
function = Link_Travel_time+transfer_penalty+Waiting_time
Searched Path Q={}
Draw Index n=1
Procedure:
Step 1: Resume Network
Step 2: Incremental n and
Step 3: Randomize edge cost (link_travel_time) by drawing the new edge cost from distribution f depends on current edge cost.
For example, new_link_travel_time:=abs(norm(link_travel_time, 5*link_travel_time)) absolute value of a normal distribution with mean = link_travel_time and standard deviation = 5* link_travel_time
Step 3: Search shortest path p_n based on randomized network
Step 3: If p_n does not exist in Q, Q=Union{Q, p_n)
Stop Criteria:
n>=N
5 Computing path attributes
Path attributes is computed based on the path composition and link attributes. Given a path consisting of a sequence of edges: {r1, r2, r3, r4….}, the following path attributes are to be computed:Total In_vehicle Travel - Time
- Total Waiting Time
- Toal Walking Time
- Total Number of Transfers
- Path-Size
5.1 Pseudo code for Computing Total In-Vehicle Time, Total Walking Time. Total Waiting Time and Total Number of Transfers
Initialization:
Total in-vehicle time = 0.0
Total Walking Time = 0.0
Total Waiting Time = 0.0
NumofTransfer=-1,
Get Path edge sequence: P={r_1, r_2, r_3,…, r_K} where K is the maximum number of edges along path
Compute index i = 1
Procedure:
Total Waiting Time = Total Waiting Time + Waiting_time(r_i)
If R_type(r_i) = “Bus” or “RTS”
Total in-vehicle time = total in-vehicle time + Link_travel_time(r_i)
NumofTransfer = NumofTransfer + 1
If R_type(r_i) = “Walk”
Total Walk Time = Total Walk Time + Link_travel_time(r_i)
Increment i
Stop Criteria:
i>K
Assign total in-vehicle time, total walking time, total waiting time and number of transfer to the path.
5.2 Pseudo code for Computing Path-Size
Computing the path-size involving checking the overlapped edges of current path with respect to other paths in the same choice set. Given a choice set C_n for OD pair n, it contains N_n paths with path i having the format of {r_1, r_2,…, r_Ki} where Ki is the maximum number of edges/route segments along path_i. The mathematical formulation of the path-size function is:
The pseudo code for computing path-size for each path i in choice set C_n is then:
Initialization:
Get Path edge sequence for each path in C_n: p_i={r_1, r_2, r_3,…, r_Ki} where Ki is the maximum number of edges along path p_i. i in the path index in choice set C_n ranging from 1 to N_n
i=1
Procedure:
Path-size=0
Path_travel_time=summation of all Link_travel_time(r) along path_i
Sub-path-size = 0 # used to store the path-size component for each edge
Sub-N = 0, # used to store the number of overlapped edge in the choice set
For each edge/route segment r along path p_i
{
Sub-path-size= Link_travel_time(r)/ Path_travel_time
For each path j in choice set C_n
{
If r exist in p_j, Sub-N = Sub-N+1
}
Sub-path-size=sub-path-size/sub-N
Path-size=path-size+sub-path-size
Assign Path-Size to path_i
}
Increment i
Stop Criteria:
i>N_n
5.3 Pseudo Code for Computing Access Walk and Egress Walk
This section presents how to compute the 1) Access walk from origin to boarding stop and 2) Egress walk from alighting stop to destination given the origin and destination from pre-day and the stop pair choice set with respect to the origin and destination.
Denote the geometry points of given OD pair is characterized by
- origin$lat for latitude of origin point
- origin$lon for longitude of origin point
- dest$lat for latitude of destination point
- dest$lon for longitude of destination point
Denote the geometry points for each stop pair {BS_i, AS_i} in stop pair choice set of given OD pair is characterized by
- BS_i$lat for latitude of boarding stop point
- BS_i$lon for longitude of boarding stop point
- AS_i$lat for latitude of alighting stop point
- AS_i$lon for longitude of alighting stop point
There are two ways to compute the access walk and egress walkthem: 1) direct distance method and 2) shortest path method
5.3.1 Method 1 – Direct Distance
For each stop pair {BS_i, AS_i} in the stop choice set of given OD:
The access walk distance dis_acess_i = geo_dist(origin$lat, origin$lon, BS_i$lat, BS_i$lon)
The egress walk distance dis_egress_i = geo_dist(AS_i$lat, AS_i$lon, dest$lat, dest$lon)
Where geo_dist is the function to compute direct distance in km between two geometry points based on the latitude and longitude of these two points.
The access walk time in minute is then time_access_i = dis_access_i/4*60
The egress walk time in minute is then time_egress_i = dis_egress_i/4*60
Store them as attributes to each stop pair in the stop pair choice set of given OD.
5.3.2 Method 2 – Shortest Path Method
Shortest path method makes use of the pedestrian network, or road network if pedestrian network is not available.
Step 1: map origin and destination to the nearest link on the road network, denote as l_origin, and l_dest
Step 2: for each stop pair {BS_i, AS_i} in the stop pair choice set for given OD:
Map BS_i and AS_i to the nearest link on the road network, denote a l_BS_i, l_AS_i.
Compute shortest path distance between l_origin and l_BS_i as access walk distance: dis_acess_i =shortest_dist(l_origin, l_BS_i)
Compute shortest path distance between l_AS_i and l_dest as the egress wlak distance: dis_egress_i = shortest_dist(l_AS_i, l_dest)
Convert distance to travel time in minutes:
time _access_i = dis_access_i/4*60
time_ egress_i = dis_egress_i/4*60
Store them as attributes to each stop pair in the stop pair choice set of given OD.
6 Route Choice Model
The pre-trip route choice model for public transportation in SimMobility is a two-level nested logit model. The first level is selection of stop pair given origin and destination and the second level is selection of stop-to-stop path given stop pair. The modelling framework is illustrated as in the Figure below.
6.1 Stop-to-Stop Path Utility given Stop pair
The utility for stop-to-stop path i in the choice set Cs given stop pair s = {boarding stop, alighting stop} will be computed based on:
where twait represents total waiting time ti vehicle, stands for total in-vehicle travel time, twalk is total walking time, ntransfer denotes the total number of transfer, β1,β2,β3,β4,β5 are the coefficients given by modeler, ln PSi is the logarithm of the path-size and represents the unobserved error component in path utility. Vi is the deterministic part of the utility which will be used to compute path selection.
In an agent-based simulation environmental, there are two ways to assign the passengers to each path in the choice set. The first approach is to follow path selection probability and assign the same portion of passengers to each path. Another way is to simulate each agent’s path selection.
6.2 Stop-to-Stop Path selection probability:
The path selection can be computed by:
6.3 Stop Pair Utility
In this benchmark model, the stop pair utility only contains three parts:
- Logsum from lower level stop-to-stop path selection
- Access walk from origin to boarding stop
- Egress walk from alighting stop to destination
TO BE ADDED
Computation for the logsum is . Please refer to the pre-day nested logit structure for detailed implementation notes.
6.4 Stop Pair selection probability:
The Stop pair path selection probability can be computed by: TO BE ADDED
Where s denotes the stop pair including initial boarding stop and final alighting stop, Cod is the stop pair choice set for OD pair od; Vs is the deterministic utility part for stop pair utility including the logsum.
6.5 Simulating passenger’s path selection
To simulate passenger’s path selection, the error component for each path in the choice that agent is facing is to be drawn from the Gumbel distribution with mean = 0, and variation = (π*π)/6. The final path with maximum utility will be chosen by that agent.
#Appendix 1 - R code for Labeling Approach
#Labeling approach on route segment network
#V3 is used for CSG test. Alone
#version1.0 convert it into a function.
# g is the igraph network
# link_travel_time stores the travel cost
#version 2.0 combined approach with K-shortest path
#update on 21ONoV13 to check the first shortest path before decising K to uses
#temp_links is also global variable
Labeling_V3<-function(origin, dest, links, stops)
{
#~ K1=1
#~ K2=6
#~ K3=12
#K=1
Big_Value<-20000000.0
# links<-temp_links[r_type!="X"]
setkeyv(links, c("start_stops", "r_type"))
links[J(origin, "Walk"), link_travel_time:=as.integer(2*Big_Value)]
#~ cat("Label 1: \n")
temp_g<-graph.data.frame(links, directed = TRUE, vertices = stops)
new_path<- get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty+E(temp_g)$waiting_time, output= "epath")
paths_hub<-list()
new_path_nodes<-Format_Path_Node(new_path,"epath", temp_g)
if(FeasiblePath(new_path, "epath", temp_g)==1 && (!(new_path_nodes %in% paths_hub)))
{
paths_hub[length(paths_hub)+1] <- new_path_nodes
#cat("Label 3 path : ", "expected cost:", compute_pathcost(new_path, "epath"), "route:",transform_path(new_path, "epath"), "\n")
}
#--------------Label 2: smallest no of transfers ------------------------------
#~ cat("Label 2: \n")
temp_links<-links[1:N_Edge]
setkeyv(temp_links, c("r_type"))
temp_links["Bus", link_travel_time:=0L]
temp_links["RTS", link_travel_time:=0L]
temp_links[,transfer_penalty:=Big_Value]
temp_links[,waiting_time:=as.integer(Big_Value)]
temp_g<-graph.data.frame(temp_links, directed = TRUE, vertices = stops)
new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty+E(temp_g)$waiting_time, output= "epath")
new_path_nodes<-Format_Path_Node(new_path,"epath", temp_g)
if(FeasiblePath(new_path, "epath", temp_g)==1 && (!(new_path_nodes %in% paths_hub)))
{
paths_hub[length(paths_hub)+1] <- new_path_nodes
#cat("Label 3 path : ", "expected cost:", compute_pathcost(new_path, "epath"), "route:",transform_path(new_path, "epath"), "\n")
}
#----------label 3: shortest walking distance---------------
#~ cat("Label 3: \n")
temp_links<-links[1:N_Edge]
setkeyv(temp_links, c("r_type"))
temp_links["Walk",transfer_penalty:=Big_Value]
temp_g<-graph.data.frame(temp_links, directed = TRUE, vertices = stops)
new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty+E(temp_g)$waiting_time, output= "epath")
new_path_nodes<-Format_Path_Node(new_path,"epath", temp_g)
if(FeasiblePath(new_path, "epath", temp_g)==1 && (!(new_path_nodes %in% paths_hub)))
{
paths_hub[length(paths_hub)+1] <- new_path_nodes
#cat("Label 3 path : ", "expected cost:", compute_pathcost(new_path, "epath"), "route:",transform_path(new_path, "epath"), "\n")
}
#~ cat("Label 4: \n")
# ------------label 4: longest MRT trips---------------------
temp_links<-links[1:N_Edge]
setkeyv(temp_links, c("r_type"))
temp_links["Walk",transfer_penalty:=Big_Value]
temp_links["Bus",transfer_penalty:=Big_Value]
temp_g<-graph.data.frame(temp_links, directed = TRUE, vertices = stops)
new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty+E(temp_g)$waiting_time, output= "epath")
new_path_nodes<-Format_Path_Node(new_path,"epath", temp_g)
if(FeasiblePath(new_path, "epath", temp_g)==1 && (!(new_path_nodes %in% paths_hub)))
{
paths_hub[length(paths_hub)+1] <- new_path_nodes
#cat("Label 3 path : ", "expected cost:", compute_pathcost(new_path, "epath"), "route:",transform_path(new_path, "epath"), "\n")
}
#-------------label 5: Minimual waiting time at transfer-------------------
#~ cat("Label 5: \n")
#~ temp_g<-g
#~ E(temp_g)[r_type!="Walk"]$link_travel_time<-0.0
#~ #new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$transfer_penalty, output= "epath"
#~ new_paths<-k_shortest_g(temp_g, origin, dest, K)
temp_links<-links[1:N_Edge]
setkeyv(temp_links, c("r_type"))
temp_links["RTS",link_travel_time:=0L]
temp_links["Bus",link_travel_time:=0L]
temp_g<-graph.data.frame(temp_links, directed = TRUE, vertices = stops)
new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty+E(temp_g)$waiting_time, output= "epath")
new_path_nodes<-Format_Path_Node(new_path,"epath", temp_g)
if(FeasiblePath(new_path, "epath", temp_g)==1 && (!(new_path_nodes %in% paths_hub)))
{
paths_hub[length(paths_hub)+1] <- new_path_nodes
#cat("Label 3 path : ", "expected cost:", compute_pathcost(new_path, "epath"), "route:",transform_path(new_path, "epath"), "\n")
}
#-----------------label 6: longest Bus trips
#~ cat("Label 6: \n")
#~ temp_g<-g
#~ E(temp_g)[r_type !="Bus"]$transfer_penalty<-Big_Value
#~ #new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty, output= "epath")
#~ new_paths<-k_shortest_g(temp_g, origin, dest, K)
temp_links<-links[1:N_Edge]
setkeyv(temp_links, c("r_type"))
temp_links["Walk",transfer_penalty:=Big_Value]
temp_links["RTS",transfer_penalty:=Big_Value]
temp_g<-graph.data.frame(temp_links, directed = TRUE, vertices = stops)
new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty+E(temp_g)$waiting_time, output= "epath")
new_path_nodes<-Format_Path_Node(new_path,"epath", temp_g)
if(FeasiblePath(new_path, "epath", temp_g)==1 && (!(new_path_nodes %in% paths_hub)))
{
paths_hub[length(paths_hub)+1] <- new_path_nodes
#cat("Label 3 path : ", "expected cost:", compute_pathcost(new_path, "epath"), "route:",transform_path(new_path, "epath"), "\n")
}
#---------Lable 7: minimal total in-vehicle travel time
#~ cat("Label 7: \n")
temp_links<-links[1:N_Edge]
setkeyv(temp_links, c("r_type"))
temp_g<-graph.data.frame(temp_links, directed = TRUE, vertices = stops)
new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty, output= "epath")
new_path_nodes<-Format_Path_Node(new_path,"epath", temp_g)
if(FeasiblePath(new_path, "epath", temp_g)==1 && (!(new_path_nodes %in% paths_hub)))
{
paths_hub[length(paths_hub)+1] <- new_path_nodes
#cat("Label 3 path : ", "expected cost:", compute_pathcost(new_path, "epath"), "route:",transform_path(new_path, "epath"), "\n")
}
#----Label 8: Minimal total travel time 2( in-vehilce+5waiting+3walking)
#~ cat("Label 8: \n")
temp_links["Walk", link_travel_time:=5L*link_travel_time]
temp_links[, waiting_time:=3L*waiting_time]
temp_g<-graph.data.frame(temp_links, directed = TRUE, vertices = stops)
new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty+E(temp_g)$waiting_time, output= "epath")
new_path_nodes<-Format_Path_Node(new_path,"epath", temp_g)
if(FeasiblePath(new_path, "epath", temp_g)==1 && (!(new_path_nodes %in% paths_hub)))
{
paths_hub[length(paths_hub)+1] <- new_path_nodes
#cat("Label 3 path : ", "expected cost:", compute_pathcost(new_path, "epath"), "route:",transform_path(new_path, "epath"), "\n")
}
#---Lable 9: Minimal total travel time 3( in-vehilce+5waiting+5walking)
#~ cat("Label 9: \n")
temp_links<-links[1:N_Edge]
setkeyv(temp_links, c("r_type"))
temp_links["Walk", link_travel_time:=5L*link_travel_time]
temp_links[, waiting_time:=5L*waiting_time]
temp_g<-graph.data.frame(temp_links, directed = TRUE, vertices = stops)
new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty+E(temp_g)$waiting_time, output= "epath")
new_path_nodes<-Format_Path_Node(new_path,"epath", temp_g)
if(FeasiblePath(new_path, "epath", temp_g)==1 && (!(new_path_nodes %in% paths_hub)))
{
paths_hub[length(paths_hub)+1] <- new_path_nodes
#cat("Label 3 path : ", "expected cost:", compute_pathcost(new_path, "epath"), "route:",transform_path(new_path, "epath"), "\n")
}
#Minimal total travel time 3( in-vehilce+5waiting+8walking)
#~ cat("Label 10: \n")
temp_links<-links[1:N_Edge]
setkeyv(temp_links, c("r_type"))
temp_links["Walk", link_travel_time:=5L*link_travel_time]
temp_links[, waiting_time:=8L*waiting_time]
temp_g<-graph.data.frame(temp_links, directed = TRUE, vertices = stops)
new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty+E(temp_g)$waiting_time, output= "epath")
new_path_nodes<-Format_Path_Node(new_path,"epath", temp_g)
if(FeasiblePath(new_path, "epath", temp_g)==1 && (!(new_path_nodes %in% paths_hub)))
{
paths_hub[length(paths_hub)+1] <- new_path_nodes
#cat("Label 3 path : ", "expected cost:", compute_pathcost(new_path, "epath"), "route:",transform_path(new_path, "epath"), "\n")
}
#label 7: lin_travel_time +2* transfer_penalty
return(paths_hub)
}
Appendix 2 - R code for Link Elimination Method
#LinkElimination_SingleOD_V2.R
#version1.0 convert it into a function.
#version 2.0 used for CSG test
# V2: links and stops as input
#V2: convert path format into nodes
#Eliminate each roadnetwork link in g and search for shortest paths
#updated on 08Oct13 to speed up for single OD path search
#Match_Table<-fread('/Users/RUI/Dropbox/Google Transit Network Data/route_seg_match_table.txt', sep=",", header = TRUE)
#Match_Table_Reg_To_Road<-fread('/home/rui/Desktop/simmobility/Google Transit Network Data/route_seg_to_Road_match_table.txt', sep = ",", header = TRUE)
# g is the igraph network
# link_travel_time stores the travel cost
#PT_sequence stores PT service lines stop sequence
LinkElimination_Road_SingleOD_V2<-function(origin, dest, links, stops)
{
#origin = "EW29"
#dest = "NS16"
#contraint 1: the first leg shouldn't be walking leg
Big_Value<-20000000.0
# links<-temp_links[r_type!="X"]
setkeyv(links, c("start_stops", "r_type"))
links[J(origin, "Walk"), link_travel_time:=as.integer(2*Big_Value)]
setkeyv(links, "edge_id")
start_g<-graph.data.frame(links, directed = TRUE, vertices = stops)
ShortestPath<- get.shortest.paths(start_g, origin, to = dest, weights = E(start_g)$link_travel_time+E(start_g)$transfer_penalty+E(start_g)$waiting_time, output= "epath")
paths_hub<-Format_Path_Node(ShortestPath,"epath", start_g)
for(i in 1:length(ShortestPath[1](/smart-fm/simmobility-prod/wiki/1)))#length(Road_Network_Path[1](/smart-fm/simmobility-prod/wiki/1))) #nrow(Match_Table))
{
#cat(i, "\\")
edgeid<-ShortestPath[1](/smart-fm/simmobility-prod/wiki/1)[i] # corresponding to start_g and links
#need to check whether it is walking leg
#~ setkeyv(Match_Table_EtoR, "edge_id") # edge_id is integer
#~ setkeyv(Match_Table_RtoE, "road_id") # road_id is text
if(edgeid<=N_Transit_Edge) # N_Transit_Edge global variable
{
#get the affected_edge_id and find the edgeid on network
road_edge_list<-lapply(str_split(Match_Table_EtoR[edgeid]$road_edge_id, "/"), as.integer)
#affected_edge_list<-lapply(str_split(Match_Table[Index_MT]$affected_edge_index, "/"), as.numeric)
# for each edge in affected_edge_list
#~ cat("update route seg", edgeid, "\n" )
# print(road_edge_list)
if(length(road_edge_list)>0)
{
for(road_edgeid in road_edge_list[1](/smart-fm/simmobility-prod/wiki/1)) #numeric
{
#get the affected_edgelist from Match_Table
temp_links<-links[1:N_Edge]
setkeyv(temp_links, c("edge_id"))
if(road_edgeid<=6669) #transit leg
{
affected_edge_list<-lapply(str_split(Match_Table_RtoE[as.character(road_edgeid)]$affected_edge_id, "/"), as.integer)
temp_links[affected_edge_list[1](/smart-fm/simmobility-prod/wiki/1), link_travel_time:=as.integer(Big_Value)]
#~ cat("------------affected_edge_list is:")
#~ print(affected_edge_list)
}
if(road_edgeid>6669) #walking_leg
{
temp_links[road_edgeid, link_travel_time:=as.integer(Big_Value)]
}
temp_g<-graph.data.frame(temp_links, directed = TRUE, vertices = stops)
new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty+E(temp_g)$waiting_time, output= "epath")
new_path_nodes<-Format_Path_Node(new_path,"epath", temp_g)
if(!(new_path_nodes %in% paths_hub))
{
paths_hub[length(paths_hub)+1] <- new_path_nodes
}
}#end for(road_edgeid in road_edge_list[1](/smart-fm/simmobility-prod/wiki/1)) #numeric
}#end if(length(road_edge_list)>0)
}
if(edgeid>N_Transit_Edge) #walk_leg
{
temp_links[road_edgeid, link_travel_time:=as.integer(Big_Value)]
temp_g<-graph.data.frame(temp_links, directed = TRUE, vertices = stops)
new_path<-get.shortest.paths(temp_g, origin, to = dest, weights = E(temp_g)$link_travel_time+E(temp_g)$transfer_penalty+E(temp_g)$waiting_time, output= "epath")
new_path_nodes<-Format_Path_Node(new_path,"epath", temp_g)
if(FeasiblePath(new_path, "epath", simulate_g)==1 && (!(new_path_nodes %in% paths_hub)))
{
paths_hub[length(paths_hub)+1] <- new_path
}
} # end if(edgeid>N_Transit_Link) #walk_leg
} # end for(i in 1:length(ShortestPath[1](/smart-fm/simmobility-prod/wiki/1)))
return(paths_hub)
}
Appendix 3 - R code for K shortest path approach
#Kshortest_edge.R
#version1.0 convert it into a function.
# g is the igraph network
# link_travel_time stores the travel cost
#update to use keys in data.table than search and update in igraph network on 21Nov13
#V2 update PT_Links to be input variable
#update to output nodes only
#add feasibility check during output
k_shortest_edge_V2<-function(origin, dest, K, links, stops)
{
#origin = "EW29"
#dest = "NS16"
#contraint 1: the first leg shouldn't be walking leg
Big_Value<-20000000.0
# links<-PT_links_InUse[r_type!="X"]
setkeyv(links, c("start_stops", "r_type"))
links[J(origin, "Walk"), link_travel_time:=as.integer(2*Big_Value)]
simulate_links<-links[1:N_links]
start_g<-graph.data.frame(simulate_links, directed = TRUE, vertices = stops)
paths_hub<-get.shortest.paths(start_g, origin, to = dest, weights = E(start_g)$link_travel_time+E(start_g)$transfer_penalty+E(start_g)$waiting_time, output= "epath")
#print(paths_hub)
N_links<-nrow(links)
simulate_g<-start_g #incase K=1
if(K>1) # added on 25Sep13 to accommandate K = 1 input
{
paths_SetB<-list()
# search for the rest
for(k in 2:K)
{
# need to address when there is only 1 single route segment in path. 25Sep13 delete -1 for compare_edges
simulate_links<-links[1:N_links]
#~ setkeyv(simulate_links, c("start_stops", "r_type"))
compare_edges<-length(paths_hub[k-1](/smart-fm/simmobility-prod/wiki/k-1))
# problem with ODs that the first shortest path is via only 1 link.
# the rest of the paths found is following the same stops pairs for the second shortest path
# updated on 03Oct, the first compare should be always eliminate the first edge of all paths
#setkeyv(simulate_links, c("edge_id"))
for(j in 1:length(paths_hub))
{
# updated on 03Oct, the first compare should be always eliminate the first edge of all paths
edgeid<-paths_hub[j](/smart-fm/simmobility-prod/wiki/j)[1]
simulate_links[edgeid, link_travel_time:=as.integer(Big_Value)]
}
for(i in 1: compare_edges) # for edgelist in path (k-1)th path in set A check against each edge in the path
{
# cat("OK2 ")
compare_list<- paths_hub[k-1](/smart-fm/simmobility-prod/wiki/k-1)[1:i]
#check whether this list happens in other shortest paths
#setkeyv(simulate_links, c("edge_id"))
for(j in 1:length(paths_hub))
{
#cat("OK3 ")
if(identical(compare_list,paths_hub[j](/smart-fm/simmobility-prod/wiki/j)[1:i]))
{
#edgeid<-paths_hub[j](/smart-fm/simmobility-prod/wiki/j)[i+1]
edgeid<-paths_hub[j](/smart-fm/simmobility-prod/wiki/j)[i] #update to [i] on 23Sep to allow change for the first link
#link_edge_id<-E(simulate_g)$edge_id
#~ vertex_id <- as.character(get.edge(simulate_g, edgeid)[1])
#~ r_type_value<-E(simulate_g)[edgeid]$r_type
#if(E(simulate_g)$r_type[edgeid]!="RTS") #added on 14Oct to give more repeativiness to MRT legs
#{
#E(simulate_g)$link_travel_time[edgeid]<-Big_Value #+E(simulate_g)$link_travel_time[edgeid]
simulate_links[edgeid, link_travel_time:=as.integer(Big_Value)]
#~ simulate_links[J(vertex_id, r_type_value), link_travel_time:=as.integer(Big_Value)]
#} #removed condition on 23Oct13
}
}
#get the start vertex for subpath searching
#start_edge_id<-paths_hub[k-1](/smart-fm/simmobility-prod/wiki/k-1)[i+1]
start_edge_id<-paths_hub[k-1](/smart-fm/simmobility-prod/wiki/k-1)[i] #update to [i] on 23Sep to allow change for the first link
#cat("start_edge_id is ", start_edge_id)
#------------Constraint 2------------
# no two walking is allowed.
#if(E(start_g)$r_service_lines[start_edge_id]=="Walk")
if(simulate_links[start_edge_id]$r_type=="Walk")
{
#get vertex node
#~ adjacent_linkid<-incident(simulate_g, get.edges(simulate_g, start_edge_id)[2], mode = "out")
#~ for(id in adjacent_linkid)
#~ {
#~ if(E(simulate_g)$r_service_lines[id]=="Walk")
#~ {
#~ #cat(id, "/")
#~ E(simulate_g)$link_travel_time[id]<-2*Big_Value#+E(simulate_g)$link_travel_time[id]
#~ }
#~ }
vertex_id <- get.edge(simulate_g, start_edge_id)[2]
start_vertex_id<-V(simulate_g)[vertex_id]$name
#~ setkeyv(simulate_links, c("start_stops", "r_type"))
if(start_vertex_id!=origin)
{
simulate_links[J(start_vertex_id, "Walk"), link_travel_time:=as.integer(Big_Value)]
}
# cat("start_vertex_id 1 is ", start_vertex_id)
}
#cat("start_vertex_id 2 is ", start_vertex_id)
simulate_g<-graph.data.frame(simulate_links, directed = TRUE, vertices = stops)
vertex_id <- get.edge(simulate_g, start_edge_id)[1]
start_vertex_id<-V(simulate_g)[vertex_id]$name
append_sub_path<- paths_hub[k-1](/smart-fm/simmobility-prod/wiki/k-1)[0:(i-1)] #start from i-1, if i = 1, this append should be empty
new_sub_path <- get.shortest.paths(simulate_g, start_vertex_id, to = dest, weights = E(simulate_g)$link_travel_time+E(simulate_g)$transfer_penalty+E(simulate_g)$waiting_time, output="epath")
new_path <- list(append(append_sub_path, unlist(new_sub_path)))
#cat("searching ", k, "th path found: ", transform_path(new_path, "epath"),"\n")
#if(length(intersect(paths_SetB, new_path))==0 & length(intersect(paths_hub, new_path))==0 ) # add new paths into paths_hub
if((!(new_path %in% paths_hub)) & (!(new_path %in% paths_SetB)))
{
#paths_hub
# cat("add new_path\n")
# print(new_path)
if(FeasiblePath(new_path, "epath", simulate_g)==1)
{
paths_SetB[length(paths_SetB)+1] <- new_path
}
# cat("new_path added\n")
}
}
#select the shortest path from SetB as kth shortest path
#cat("OK5 ")
if(length(paths_SetB)==0)
{
#cat("at ",k, " breaked \n")
break;
}
#
for(l in 1: length(paths_SetB))
{
#cat("OK6 ")
if(l ==1)
{
min_path<- compute_pathcost_V2(paths_SetB[1], "epath", simulate_g)#sum(E(g)$link_travel_time[unlist(paths_SetB[1])])
min_path_id <- 1
}
current<-compute_pathcost_V2(paths_SetB[l], "epath", simulate_g)#sum(E(g)$link_travel_time[unlist(paths_SetB[l])])
if(current<min_path)
{
min_path<- current
min_path_id <- l
}
}
#cat("OK7 ")
#update new_path and put into paths_setB
paths_hub[k]<-paths_SetB[min_path_id]
#print(paths_hub[k])
#paths_hub[k] <- Format_Path_Node(paths_SetB[min_path_id],"epath", simulate_g)
paths_SetB[min_path_id]<-NULL
}
}
#cat("--------------path no: node, service lines on route_seg, node, service lines on route_seg, node---------------- \n")
#for( i in 1: length(paths_hub))
#{
# cat("path ", i, ": ", "expected cost:", compute_pathcost(paths_hub[i], "epath"), "route:",transform_path(paths_hub[i], "epath"), "\n")
#}
#paths_hub<-Format_Path_Node(new_path,"epath", start_g)
if(FeasiblePath(paths_hub[1], "epath", start_g)==0)
{
paths_hub[1]<-NULL
}
for( i in 1: length(paths_hub))
{
paths_hub[i]<-Format_Path_Node(paths_hub[i],"epath", simulate_g)
#check feasibility of paths_hub[i]
}
return(paths_hub)
}
Appendix 4 - R code for simulation approach
# this is simulation approach
# Simulation approach is going to add variance to each link in the network and re-search for shortest path
# propose to use lognormal distribution. mean? variance?
# network i iGraph format as Input_NW
# -----lognormal--------------
# dlnorm(x, meanlog = 0, sdlog = 1, log = FALSE)
#plnorm(q, meanlog = 0, sdlog = 1, lower.tail = TRUE, log.p = FALSE)
#qlnorm(p, meanlog = 0, sdlog = 1, lower.tail = TRUE, log.p = FALSE)
#rlnorm(n, meanlog = 0, sdlog = 1)
#N is the total number of draws
Simulation<-function(origin,dest,N, links, stops){
Big_Value<-20000000.0 #increase this value to be larger then the big_value in labeling_V2/V1 on 23Oct to avoid effect concellation
Big_Value<-20000000.0
setkeyv(links, c("start_stops", "r_type"))
links[J(origin, "Walk"), link_travel_time:=as.integer(Big_Value)]
start_g<-graph.data.frame(links, directed = TRUE, vertices = stops)
ShortestPath<- get.shortest.paths(start_g, origin, to = dest, weights = E(start_g)$link_travel_time+E(start_g)$transfer_penalty+E(start_g)$waiting_time, output= "epath")
paths_hub<-Format_Path_Node(ShortestPath,"epath", start_g)
N_links<-nrow(links)
if(N>1) #
{
#~ counter = 2
for(n in 2:N)
{
simulate_links<-links[1:N_links]
simulate_links$link_travel_time<-abs(rnorm(N_links, mean= simulate_links$link_travel_time, 5*simulate_links$link_travel_time))
simulate_g<-graph.data.frame(simulate_links, directed = TRUE, vertices = stops)
#~
#simulate_g<-start_g
#set.seed(n)
#E(simulate_g)$link_travel_time<-E(simulate_g)$link_travel_time + dlnorm(1, meanlog = 0, sdlog = 10/E(simulate_g)$link_travel_time)
#E(simulate_g)$transfer_penalty<-Big_Value
new_path<-get.shortest.paths(simulate_g, origin, to = dest, weights = E(simulate_g)$link_travel_time+E(simulate_g)$transfer_penalty+E(simulate_g)$waiting_time, output= "epath")
new_path_nodes<-Format_Path_Node(new_path,"epath", simulate_g)
if(FeasiblePath(new_path, "epath", simulate_g)==1 && (!(new_path_nodes %in% paths_hub)))
{
paths_hub[length(paths_hub)+1] <- new_path_nodes
} # end if(FeasiblePath(new_path, "epath", simulate_g)==1 && (!(new_path_nodes %in% paths_hub)))
}
}
#~ #cat("--------------path no: node, service lines on route_seg, node, service lines on route_seg, node---------------- \n")
#~ for( i in 1: length(paths_hub))
#~ {
#~ cat("path ", i, ": ", "expected cost:", compute_pathcost(paths_hub[i], "epath"), "route:",transform_path(paths_hub[i], "epath"), "\n")
#~ }
return(paths_hub)
}