Optimization of Public Transport Services to Minimize Passengers’ Waiting Times and Maximize Vehicles’ Occupancy Ratios

Abstract

Determining the best timetable for vehicles in a public transportation (PT) network is a complex problem, especially because it is just necessary to consider the requirements and satisfaction of passengers as the requirements of transportation companies. In this paper, a model of the PT timetabling problem which takes into consideration the passenger waiting time (PWT) at a station and the vehicle occupancy ratio (VOR) is proposed. The solution aims to minimize PWT and maximize VOR. Due to the large search space of the problem, we use a multiobjective particle swarm optimization (MOPSO) algorithm to arrive at the solution of the problem. The results of the proposed method are compared with similar results from the existing literature.

1. Introduction

In modern transportation systems, the greatest challenge is to minimize, in general terms, energy consumption, and maximize economic, technological and social goals. The problem of optimizing and finding the best timetable for public transportation (PT) vehicles has been known for years [1,2,3]. Recent research has provided more efficient algorithms which have achieved better results by modifying known mathematical models and modifying or combining various known algorithms [4,5,6]. Planning of PT is a highly complex task which is usually analyzed via two different aspects: minimization of the passenger waiting time (PWT) at the station and optimization of the number and/or sizes of vehicles. The train timetabling problem (TTP) has recently been studied, and the main problem in this field is to determine a periodic or non-periodic timetable, which satisfies the capacities of vehicles and limits of operations [7,8,9,10]. Some of the greatest challenges in waiting time (WT) minimization models are to find the optimal number of vehicles and to find the optimal route and minimize travel time [11]. The goal of every PT service should be to attract more people to use it by reducing the use of private cars, which is directly related to reducing traffic congestion, decreasing the number of car accidents and reducing pollution. The use of PT services by passengers depends on three elements; namely, travel time (walking, waiting and riding times), fare (ticket and other related services’ costs) and convenience (a comfortable walk, waiting under a shelter, having a seat on the vehicle, air-conditioning on the vehicle, etc.). However, for the PT operator to maximize profits, operational costs need to be minimized. The PT operator’s requirements can be achieved by designing an efficient network (the more transfers the network has the more efficient it is), adopting a quality timetable, efficiently schedule the vehicles and maximize the vehicle occupancy ratio (VOR). This paper proposes a model of the PT timetabling problem whose solution improves PT operations planning in terms of timetabling and vehicle scheduling; that is, changes to the departure times and assignment of the PT vehicles, so as to reduce total PWT and increase VOR (thereby minimizing operational costs of the vehicles for the PT operator). Due to the complexity of the problem, a multiobjective particle swarm optimization (MOPSO) algorithm is used in order to minimize PWT while maximizing VOR. The paper is structured as follows: Section 2 extensively describes the state-of-the-art in the timetabling problem. Section 3 elaborates on the problem statement. Section 4 presents the proposed method, and Section 5 contains two numerical examples. Section 6 presents the analysis and discussion of the results. Section 7 contains the concluding remarks.

2. State-of-the-Art

In order to increase the productivity and efficiency of transport services on the one hand, and customer satisfaction on the other, four main analytical methods for determining the number of vehicles needed during the relevant period are presented in [1,3,12]. The first two methods in these papers ensure the average maximum daily occupancy of the vehicles during a given period. The other two methods elaborate on the capacities of the vehicles, which will never be exceeded with an additional constraint on the part of the route which is loaded more than the required availability. Numerical calculation of the average PWT at a station with a limited capacity of the intercity transport vehicles is described in [13]. It is concluded that for a more accurate representation of PWT at the station, the reliability of supply services, passenger behavior and characteristics of the transportation system should be considered. The passenger transport problem is always observed from the viewpoint of costs and improvement of business transport companies. On the contrary, passenger satisfaction is rarely taken into account, especially when creating the timetable [14]. Passenger satisfaction is directly related to the total travel time; the shorter the travel time the more satisfied the passenger is and vice-versa. Excluding the travel speed of the vehicle, the total travel time can be reduced by reducing PWT and the number of passenger transfers. Hence, in this paper, we use the term passenger satisfaction to refer primarily to PWT. A multicriteria approach to timetabling problems, with an emphasis on minimizing PWT at the station, is proposed in [15,16]. The authors in these papers analyze two criteria; namely, empty seat penalty (empty seat kilometers or empty seat hours) and approximate PWT at the station. The problem is solved using a multiobjective label-correcting algorithm that results in 43% saving of PWT with an acceptable load of the vehicle. In order to ensure fast and energy efficient urban rail transport, a nonlinear problem of minimizing the total travel time and energy consumption is presented in [17]. The model reduces the cost (number of trains and energy consumption) and improves passenger satisfaction (reduced PWT and the number of transfers, thereby reducing the total travel time). Optimization of energy consumption and PWT can also be found in [18]. The authors in this paper propose a bi-objective timetable optimization model to minimize PWT and energy consumption. The genetic algorithm (GA) is used for generating a solution with reduced total energy consumption and total passenger waiting time in comparison to the real timetable. The timetable synchronization problem (TSP) refers to the problem of waiting time during a transfer, which is usually solved using a branch and bound (B&B) method. In order to speed up the execution time, the optimization based heuristic method (OHM) is developed and compared with B&B in [19]. It is concluded that OHM is much faster and optimizes the problem more efficiently. The scheduling problem is usually based on maximization of the number of synchronized vehicles arriving at the transfer station or minimization of the total waiting time at the station. For solving the latter problem, a genetic algorithm with local search is used in [20]. The model is applied to a small bus network and the cost of the waiting time is reduced by 9.5%. An optimization model for synchronizing a timetable is proposed in [21]; i.e., for minimizing the maximum PWT during the transfer and reducing the worst time of the transfer. Mathematical models with PWT as the objective function at the transfer station include a set of mixed integer programming (MIP) models [22]. The model can be solved using a conventional MIP solver such as CPLEX solver (B&B) if there are fewer than 50 lines and by using genetic algorithms for a greater network. Additionally, a solution of the timetable synchronization problem is proposed in [23]. The authors develop a multicriteria optimization model which takes the vehicle scheduling and passenger demand assignment into account. In order to solve the problem to obtain a set of Pareto-efficient solutions, a novel deficit function (DF)-based sequential search method is proposed. Minimization of the average transfer time in the periodic scheduling of trains (PRTS—periodic railway timetable scheduling problem) is solved using an improved differential evolution (DE) algorithm with dual population [24]. A comparison of the presented model against the B&B method and greedy-based heuristic algorithm for using the PRTS simulation algorithm to solve the problem of schedules shows that the given model provides a better indicator for PRTS problem and numerical indicators of optimization functions. The problem of minimizing PWT at the station and the cost of vehicle occupancy is solved using the genetic algorithm in [25]. The models developed for solving the problem of minimizing PWT at the station are defined with preserving the flow of passengers waiting at the station, as shown in [26,27]. The problem of PWT at the station occurs when there is a delay in the timetable. One possible solution is to include additional time delays according to the probability theory; i.e., the objective function includes an exponential distribution of delay, $Exp\left(\right)$ with the expected delay, as shown in [28,29]. In [30], the problem of minimizing train delays while maximizing the total satisfaction during a traffic jam (for example, a vehicular crash or similar) is solved using a heuristic algorithm. Timetable optimization is based on the optimal departure time of vehicles from the station for each line of vehicles in order to reduce PWT. A model which minimizes PWT and distinguishes a direct vehicle transfer from walking from one station to another is solved using a heuristic algorithm, the same as in [31]. If the headway is reduced, the average PWT can be reduced as well. Minimization of the sum of headways results in minimization of the average PWT. The average PWT is equal to the ratio of the sum of headway squares and the sum of headways during which a traveler arrives. A numerical method for solving this issue is described in [32].

In order to rationalize departure intervals as much as possible, make the bus journey quicker and minimize PWT, an optimization bus schedule model is proposed in [33]. The authors in this paper consider vehicle overtaking, the limit of the vehicle’s capacity and the uncertainty of passenger choice for a bus type (traditional bus or rapid bus). They propose two methods for the solution: a hybrid method of traditional PSO (HPSO) and a combination of GA and PSO named GAPSO.

Table 1. Literature review.

Paper	Year	Objective	Constraints	Solution Method	Case
[28]	2006	Minimize the waiting time	Running time, departure time, buffer time, spacing time	Linear programming	real
[19]	2008	Minimize the interchange waiting times	Running time, dwell times, trip times, headways, turnaround	B&B, heuristic	real
[22]	2010	Minimize the waiting time	Headway bounds, extra dwell times	Genetic algorithm	real
[15]	2012	Minimize the expected passenger waiting time and Minimize the discrepancy from a desired occupancy level on the vehicles	Headway bounds, vehicle capacity	A multi-objective label-correcting algorithm	real
[25]	2012	Minimize the waiting time	Headway, in-train passengers, passenger demands, fleet-size, number of boarded passengers	Genetic algorithm	real
[24]	2013	Minimize the waiting time	Traveling time	Evolutionary algorithm	real
[31]	2014	Minimize the waiting time	Headway, departure time	Tabu Search algorithm	test
[17]	2015	Minimize the total travel time of all passengers and the energy consumption of the trains using a weighted sum strategy	Headway, train capacity	Evolutionary algorithms	test
[21]	2015	Minimize the maximal passenger waiting time	Headway, departure time, running time	Genetic algorithm	real
[32]	2015	Minimize the waiting time	Headway	Analytical	test
[33]	2017	Minimize the waiting time	Headway	Improvements of the Genetic Algorithm and Particle Swarm Optimization	test
[23]	2018	Vehicle scheduling problem with the transit assignment	Headway, the number of vehicle departures, fleet size	Heuristic	test
[18]	2019	Minimize the PWT and energy consumption	Headway	Genetic algorithm	real
This paper		Minimize the waiting time and maximize vehicles’ occupancy	Headway, passenger demands, number of boarded passengers	Particle swarm optimization	test

summarizes a literature review and shows the details of the studies presented in this section.

3. Problem Statement

In order to further elaborate on the effects of departure times and the assignments of PT vehicles on PWT and VOR, a representative network with stations A, B, C, D and E is considered [23] (see Figure 1). The PT network has two terminals (a and b), two routes ($r_{a\stackrel{}{\to }b}$, $r_{b\stackrel{}{\to }a}$ and one transfer stop (node C). An estimated origin-destination (OD) demand matrix is shown in

Table 2. Simple example origin destination matrix [23].

From\To	A	B	C	D	E
A	0	100	50	70	80
B	0	0	0	0	50
C	50	20	0	80	60
D	60	0	0	0	0
E	80	100	20	60	0

.

The maximum load on route $r_{a\stackrel{}{\to }b}$ is the load on route segment C-B (360 passengers), and the maximum load on route $r_{b\stackrel{}{\to }a}$ is the load on route segment C-D (340 passengers). Assuming that a vehicle with a suitable capacity is used for service ACBE1, all 360 passengers on route segment C-B will be satisfied; i.e., transferred from their respective origins to destinations. For a case in which the vehicle does not have a suitable capacity, some passengers will be left at one of the stations, depending on the vehicle capacity. These passengers may choose to wait for the next service or use another type of transport (e.g., taxi). This decision is based upon the time needed to wait for the next service. In order to perform some example calculations, let us assume that the desired occupancy ($d_o$) for both routes is set to 70 passengers and that we have a given set of the number of departures per hour (four for route $r_{a\stackrel{}{\to }b}$ and four for route $r_{b\stackrel{}{\to }a}$). After service ACBE1 leaves Station A, 230 passengers are left behind and have to wait for the next service. Assuming that all 50 passengers for Station C get onboard at Station A, after service ACBE1 leaves Station C, 60 passengers are left behind. Considering the number of passengers left in Station A by the previous service (ACBE1), after service ACBE2 leaves Station A, 460 passengers are left behind. Assuming that all 50 passengers for Station C get onboard at Station A, after service ACBE2 leaves Station C, 120 passengers are then left behind (this number includes those passengers left behind from the previous service). The vehicle occupancy ratio for services ACBE1 and ACBE2 is the same: [1 1 1]. Each element of the occupancy ratio array is given by the ratio of the number of passengers in the vehicle to the vehicle capacity and is defined for each OD pair of the service—in this case AC, CB, BE. The amount of PWT for a given station is defined as the product of the number of passengers left behind by the previous service and the time needed for the current service or vehicle to arrive (plus an additional constant term which takes into consideration the number of passengers arriving at the station in the meantime—see Equation (11)). Hence, if the time between consecutive services is 15 min, i.e., service ACBE2 leaves 15 min after ACBE1, the amount of PWT at Station A for service ACBE2 is 3450 $passengers·min$, while that of Station C is 900 $passengers·min$. Assuming that the headway between services is too long, it can be assumed that the remaining passengers will find other means of transportation; this is suitable neither for the service operator nor for the passengers. A cost effective solution for the service operator would be if route ACBE were to be such that the same vehicle and crew can return and perform both services within an hour. Otherwise, if the operator sends more vehicles per time period, additional costs are incurred (more vehicles and crew members are needed). Supposing the operator has passenger cars or coaches that can be connected (assuming this is a tram or rail line), then, in this example, if the capacity is 350 ($5·70$), no passengers are left behind at Station A, which is good. However, the vehicle occupancy ratio for both ACBE1 and ACBE2 is now $\left[0.857110.54286\right]$. The numbers of passengers left at the stations are $\left[0100\right]$ and $\left[0200\right]$ respectively for the two consecutive services. For these reasons, it is necessary to make optimal decisions in order to satisfy the operators business interests and passenger needs at the same time.

4. Proposed Model

In order to describe the proposed model, an overview of the notation and variables used is provided in

Table 3. Notation.

GENERAL
OD	Origin-destination
$O{D}_s$	Set of OD pairs for a given service s
$S_{OD}$	Set of services for a given OD pair
$M=\{1,2,\dots ,m\}$	Set of vehicles
$N=\{1,2,\dots ,n\}$	Set of stations
L	Set of lines
$v_{s,c}$	Vehicle v with capacity c, for service s
INDEXES
$s\in S_{OD}$	Service s from the set of services
$i\in N$	Station i from the set of stations
$v\in M$	Vehicle v from the set of vehicles
$l\in L$	Line l from the set of lines
VARIABLES
$t_{s,i}^d$	Departure time of the vehicle at station i for service s
$t_{s,i}^a$	Arrival time of the vehicle at station i for service s
$d_o$	Desired occupancy
$dwel{l}_{s,i}$	Dwelling time at station i for service s
$ru{n}_{s,i}$	Running time – traveling time between adjacent
	stations i and $i+1$ for service s
$ru{n}_s$	Sum of running time between all adjacent
	stations for service s
$H_{i,s}$	Headway – difference between departure times at station i
	for consecutive services s and $s+1$
$T_{s,i}$	Difference between departure times of vehicle v
	between adjacent stations i and $i+1$ for a given service s
$T_i$	Time horizon for i-th station
$P(v_{s,c},i)$	Total number of passengers in vehicle v in station i
$P_{ij}$	Number of passengers entering the vehicle in station i
	heading for destination j ($ij$ OD pair)
$P(v_{s,c},kj)$	Number of passengers who entered the vehicle
	in station k with traveling to j
$P^{out}(v_{s,c},i)$	Number of passengers exiting vehicle v of service s at station i
$P^{in}(v_{s,c},i)$	Number of passengers entering vehicle v of service s at station i
$P^{curr}(v_{s,c},i)$	Number of passengers in vehicle v of service s arriving at station i
$P^{stay}(v_{s,c},i)$	Number of passengers remaining at station i
	after vehicle v of service s leaves the station
$w_i$	Average waiting time at station i
$k_i$	Average number of passengers per time
$Z_{i,s}$	Amount of PWT at station i for service s
$\eta$	Total vehicles’ occupancy ratio for all services for a given line
${\tau }_{v,s}$	Average vehicle occupancy ratio for service s

.

The OD demand formulation of the number of passengers is given as input data in order to calculate PWT for a given station. According to Figure 1, the OD pairs are AB, AC, AD, AE, CD, CA, CB, CE, BE, EC, EB, ED, EA and DA; while the possible lines are AC, ACB, ACBE, ACD, etc. Each line $\left(l\right)$ has multiple services per day; for example, line ACBE has services ACBE1, ACBE2, ACBE3, etc. Hence, a set of services is defined for each OD pair; i.e., $S_{ACBE}=ACBE1,ACBE2,ACBE3,\dots$ Each service s depends on the:

• Departure time of the vehicle from station i $\left(t_{s,i}^d\right)$;
• Vehicle assignment—$v_{s,c}=$ the vehicle v of service s with given capacity c.

The set of OD pairs for each service s is denoted as $O{D}_s$, and the set of services for each OD pair is denoted as $S_{OD}$.

Let us suppose that the parameters are as follows. There are six OD pairs for one line ACBE: AC, AB, AE, CB, CE and BE (i.e., $O{D}_{\left\{ACBE\right\}}=\{AC,AB,AE,CB,CE,BE\}$). The number of passengers entering service s at station i $\left(P^{in}(v_{s,c},i)\right)$ is lower than or equal to the sum of the number of passengers for all OD pairs in station i at departure time t in service s:

$${\displaystyle P^{in}(v_{s,c},i)\le \sum _{j=i+1}^n{P}_{ij}}$$

(1)

where $P_{ij}$ is the number of passengers who arrive at station i and are traveling to station j (input data from OD matrix). An inequality sign in Equation (1) signifies that is not necessary that all passengers at station i board vehicle v.

4.1. Assumptions, Variables, Parameters Which Are Time Dependent

First, a set of vehicles with different capacities is considered. The set of vehicles is marked with $M=1,2,\dots ,m$ (v-th vehicle is marked with index v). This is important in order to track the number of used vehicles in the network at a certain time. The time between two consecutive services in station i with the same OD pair is defined as:

$$H_{i,s}=t_{s,i}^d-t_{(s-1),i}^d,\forall s\in S_{OD}$$

(2)

which is the objective of timetable generation and is a dependent variable of the model.

When the vehicle arrives at the final station of the network for one service, this vehicle is free to be used for another service. The traveling time between adjacent stations is defined as running time $\left(ru{n}_{s,i}\right)$ for each vehicle v:

$$ru{n}_{s,i}=t_{s,(i+1)}^a-t_{s,i}^d.$$

(3)

For each vehicle v at each station i of the service s, the dwelling time $\left(dwel{l}_{s,i}\right)$ is defined:

$$dwel{l}_{s,i}\le t_{s,i}^d-t_{s,i}^a.$$

(4)

The time for boarding or alighting is considered as a constant for now. The difference between departure times of vehicle v between consecutive stations (running time + dwelling time) for a given service s (see Figure 2) is defined by:

$$T_{s,i}=t_{s,i+1}^d-t_{s,i}^d,\forall s\in S_{OD}.$$

(5)

4.2. Assumptions, Variables, Parameters and Sets—Passengers, PWT and VOR

For each OD pair, we have passenger variables that depend on departure times for stations $\left(t_{s,i}^d\right)$. The number of passengers in i-th station for j-th OD-pair (e.g., AB) for all stations for one line (at departure time $t_{s,i}^d$), for a given vehicle $v_{s,c}$ is shown in Figure 2. The total number of passengers in vehicle v in station i is defined by:

$$P(v_{s,c},i)=P^{curr}(v_{s,c},i)+P^{in}(v_{s,c},i)-P^{out}(v_{s,c},i)$$

(6)

where $P^{out}(v_{s,c},i)$ is the number of passengers leaving service s at station i and $P^{in}(v_{s,c},i)$ is the number of passengers entering service s at station i. The following assumptions are made for the number of passengers entering service s at station i: (a) The maximum number of passengers that can board the service depends on the currently available capacity of the service; i.e., the current number of free seats. (b) Passengers board the vehicle at station i according to the direct proportional distribution determined by the current number of passengers for each OD pair of station i (station i being the origin).

Those assumptions are necessary in order to simplify the calculation of Equation (6); i.e., $P^{out}$ for a given station will always be the maximum possible for that station.

The number of passengers in vehicle v before station i (entered in station k) and with destination j (curr = current) is:

$$P^{curr}(v_{s,c},i)=\sum _{k=1}^{i-1}\sum _{j=i+1}^{n}P(v_{s,c},kj)$$

(7)

$$P(v_{s,c},i)=\sum _{k=1}^{i-1}\sum _{j=i+1}^{n}P(v_{s,c},kj)+P^{in}(v_{s,c},i)-P^{out}(v_{s,c},i);$$

(8)

i.e., passengers entering before station i and exiting after station i (station i excluded).

The average number of passengers at a station between consecutive services.

$$k_i={\displaystyle \frac{{\sum }_{s\in S_{OD}}P(v_{s,c},i)}{T_i}}$$

(9)

where $T_i$ is the time horizon for station i. The average waiting time at station i is given by [23]

$$w_i={\displaystyle \frac{E\left[H_i\right]}{2}}\left[1+{\displaystyle \frac{Var\left[H_i\right]}{E^2\left[H_i\right]}}\right]$$

(10)

where $E\left[H_i\right]$ and $Var\left[H_i\right]$ are, respectively, the mean and variance of headway time between vehicles at station i. The amount of PWT in station i for every service (s) is given by

$$Z_{i,s}=P^{stay}(v_{s,c},i)·H_{i,s}+k_i{H}_{i,s}{w}_i,$$

(11)

where $P^{stay}(v_{s,c},i)$ is the number of passengers left at station i because they could not enter service $s-1$. The amount of PWT (Equation (11)) at station i is the sum of the amount of waiting time at station i until the next service arrives at the station $(P^{stay}(v_{s,c},i)·H_{i,s}))$ and the average amount of waiting time at station i given by $k_i{H}_{i,s}{w}_i$ (for all OD pairs). It should be noted that this constant term is a measure of passengers arriving at the station between services with the average waiting time $\left(w_i\right)$ for these passengers. This is based on the assumptions that (a) passengers randomly arrive at the station and (b) the arrivals of vehicles are uneven but occur in predetermined time intervals.

In the multiobjective optimization algorithm implemented in this paper, average normalized values of PWT and VOR are used. The average normalized amount of PWT is defined by: ):

$$Z_{av,norm}={\displaystyle \frac{{\sum }_{i,s}{Z}_{i,s}^{norm}}{n·m}}$$

(12)

where the normalized amount of PWT is defined by $Z_{i,s}^{norm}={\displaystyle \frac{Z_{i,s}-Z_{min}}{Z_{max}-Z_{min}}}$ for a given station i and service s, with $Z_{max}$ being the amount of PWT for the maximal difference between departure times of consecutive services for a given station using vehicles of minimal capacity, and $Z_{min}$ is the amount of PWT for the minimal differences between departure times of consecutive services for a given station using vehicles of maximal capacity. n is the number of stations and m is the total number of vehicles used (i.e., number of services). Equation (12) represents the objective function for minimizing the amount of PWT.

The vehicle occupancy ratio for service s is defined by:

$${\tau }_{v,s}={\displaystyle \frac{{\sum }_{i=1}^n{P}^{curr}(v_{s,c},i)·ru{n}_{s,i}}{n·c_v·ru{n}_s}}$$

(13)

where $ru{n}_s={\sum }_{i=1}^n{r}_{s,i}$ is the sum of running time between all adjacent stations for service s.

The total vehicle occupancy ratio (VOR) for all services across all stations is given by:

$$\eta =\sum _{s\in S}{\tau }_{v,s}.$$

(14)

The average normalized value of $\eta$, ${\eta }_{av,norm}$ is defined by:

$${\eta }_{av,norm}={\displaystyle \frac{\eta }{n·m}}$$

(15)

Equation (15) represents the objective function for maximizing VOR. The profit margin of a transport company increases as the vehicles’ occupancy ratio increases.

4.3. Objective Function

The optimization problem aims at minimizing the amount of PWT while maximizing VOR. The variable PWT is used in the optimization model in order to satisfy users of PT. On the another hand, VOR (Equation (15)), used in the optimization model, is based on vehicles’ occupancy for all services in time horizon. The variables included in the model are discrete so the model can be solved using combinatorial optimization methods.

The objective function of the model is defined as:

$$\underset{x}{min}F\left(x\right)=(f_1\left(x\right),f_2\left(x\right))$$

(16)

where x is decision variable in the solution space of dimension $2|S_{OD}|$. $|S_{OD}|$ represents the number of services in $S_{OD}$. The first $|S_{OD}|$ elements in x represent the departure times of the services from the first station on the line, while the second $|S_{OD}|$ elements in x represent the corresponding vehicle capacities of the services. $f_1\left(x\right)$ corresponds to the inverse average normalized value of $\eta$, i.e., $(1-{\eta }_{av,norm})$ given by (15), and $f_2\left(x\right)$ corresponds to the average normalized amount of PWT, $Z_{av,norm}$, given by (12). The objective function (16) minimizes two goals and it is solved using multiobjective optimization. It minimizes $f_1\left(x\right)$ and $f_2\left(x\right)$, which correspond to maximizing VOR and minimizing the amount of PWT. The objective function of the model has the following constraints:

$$\begin{array}{ccc}\hfill P(v_{s,c},i)\le c_v,& \hfill & \hfill \forall v\in \mathcal{M},\forall i\in N\end{array}$$

(17a)

$$\begin{array}{ccc}\hfill t_{s-1,i}^d<t_{s,i}^d,& \hfill & \hfill \forall s\in S_{OD}\end{array}$$

(17b)

$$\begin{array}{ccc}\hfill T_{s,i}=t_{s+1,i}^d-t_{s,i}^d,& \hfill & \hfill \forall i\in N,\forall s\in S_{OD}\end{array}$$

(17c)

$$\begin{array}{ccc}\hfill H_{i,s}=t_{s,i}^d-t_{s-1,i}^d,& \hfill & \hfill \forall i\in N,\forall s\in S_{OD}\end{array}$$

(17d)

$$\begin{array}{ccc}\hfill k_i\le {\displaystyle \frac{c_v}{{\displaystyle \sum _{s\in S}{H}_{i,s}}}}& \hfill & \hfill \forall i\in N,\forall s\in S_{OD}\end{array}$$

(17e)

In order to simplify the model and future calculations, it is assumed that the number of passengers in j-th vehicle cannot be more than $c_v$, $\forall v\in M$, as shown in constraint (17a). Constraint (17b) implies that the departure time of vehicle v in service s has to be after the previous departure time of the same OD pair served. Equation (17c) expresses the difference between departure times of vehicle v between consecutive stations in one service s. Equation (17d) expresses the difference of the departure time between two consecutive services in station i with the same OD pair. Constraint (17e) expresses that the average arrival rate for the OD pair at station i is less than or equal to the average maximum capacity rate.

5. Results

Due to the complexity of the objective function (16) and the large search space, we propose to use a heuristic optimization algorithm to determine a suitable solution. Such algorithms have shown to be suitable in such optimization problems, especially those involving large search space [34,35,36,37]. The efficiency and suitability of the various heuristic optimization algorithms in solving a whole range of complex problems have been receiving a lot of attention from academia for many years. Most of the available heuristic optimization algorithms mainly fall into two categories—swarm intelligence algorithms and evolutionary algorithms. The main representatives of these categories are particle swarm optimization (PSO) and the genetic algorithm (GA), respectively. In this paper, the multiobjective particle swarm optimization (MOPSO) algorithm, proposed in [38], was used. The algorithm in [38] extends the standard PSO algorithm to solve multiobjective optimization problems by utilizing an additional repository of particles to help the main set of particles in their search. The exploratory capabilities of the algorithm are also enhanced by a specific mutation operator that is incorporated. A preliminary comparison by the authors of this paper, of an implementation of the aforementioned algorithm, was made with a brute force search for the Pareto-optimal set solution of a model of a much simpler PT timetabling problem, and the same results were obtained. As a result, the implemented version was deemed suitable enough. A comparison of the MOPSO algorithm implemented in this paper with other PSO algorithms for multiobjective optimization, such as that proposed in [39], will be considered in future research, especially with respect to the accuracy and convergence rate in solving the proposed objective function. Additionally, a comparison of other heuristic optimization algorithms regarding their suitability in solving the proposed objective function will also be considered.

The proposed optimal solution from the Pareto-optimal set, obtained using the implemented MOPSO algorithm, was determined based on the multicriteria decision making method. The technique for order of preference by similarity to ideal solution (TOPSIS), as proposed in [40], was used for solving traffic problems. In this section, two experimental problems based on input data and assumptions given in [23] are provided and analyzed. In both experiments, the results obtained are compared using the proposed method and assumptions and the results from [23].

5.1. Experiment 1

The first experiment consists of a simple passenger transportation line (Figure 1) involving three sets of a number of departures (q = 4, q = 5 and q = 6) for route $r_{a\stackrel{}{\to }b}$ with four stations and two sets of a number departures (q = 4 and q = 5) for route $r_{b\stackrel{}{\to }a}$. The input data during the given time period (7 a.m.–8 a.m.) consist of the average travel times, and the constant time needed for alighting and boarding is set to 0.5 minutes. An estimated OD demand matrix is presented in Table 2. The running times between each station for each service are presented in Figure 1. Other details of the experimental setup, which define the search space, are as follows:

• Desired occupancy of each vehicle: 70;
• Possible departure time (in minutes) at the first terminal is defined by the time intervals for each service respectively:

  -

  $[0,15],[16,30],[31,45],[46,60]$

  -

  $[0,12],[13,24],[25,36],[37,48],[49,60]$

  -

  $[0,10],[11,20],[21,30],[31,40],[41,50],[51,60]$.

Table 4. Comparison of the results obtained in Experiment 1 using the proposed method and those obtained in the literature [23].

${\mathit{d}}_{\mathit{o}}=70$			Dep. Time	${\mathit{Z}}_{\mathit{av},\mathit{norm}}$	Number of Passengers Left at the Station			Waiting Time			Amount of PWT (Equation (11))
$r_{a\stackrel{}{\to }b}$	q = 4	Proposed method	[7:15, 7:28, 7:38, 7:46]	0.142	230	152	8	13	12	12	4386.45	2813.87	336.74
					460	304	16	10	9	9	5674.19	3684.52	339.03
					690	456	24	8	7	7	6379.35	4163.61	335.23
					920	608	32	90	89	89	92,467.74	60,520.65	4491.29
		Results according to [23]	[7:15, 7:30, 7:45, 8:00]	0.161	230	152	8	15	14	14	5700	3630	495
					460	304	16	15	14	14	9150	5910	615
					690	456	24	15	14	14	12,600	8190	735
					920	608	32	90	89	89	96,300	62,820	5130
	q = 5	Proposed method	[7:12, 7:24, 7:35, 7:43, 7:49]	0.126	230	152	8	12	11	11	4239.73	2711.84	342.62
					460	304	16	11	10	10	6416.42	4157.85	402.07
					690	456	24	8	7	7	6506.49	4239.89	356.41
					920	608	32	6	5	5	6259.86	4091.92	315.31
					1150	760	40	72	71	71	91,678.38	60,047.03	4359.73
		Results according to [23]	[7:12, 7:24, 7:36, 7:48, 8:00]	0.141	230	152	8	12	11	11	4560	2904	396
					460	304	16	12	11	11	7320	4728	492
					690	456	24	12	11	11	10,080	6552	588
					920	608	32	12	11	11	12,840	8376	684
					1150	760	40	72	71	71	93,600	61,200	4680
	q = 6	Proposed method	[7:10, 7:20, 7:30, 7:40, 7:47, 7:51]	0.114	230	152	8	10	9	9	3635.37	2321.22	302.56
					460	304	16	10	9	9	5935.37	3841.22	382.56
					690	456	24	10	9	9	8235.37	5361.22	462.56
					920	608	32	7	6	6	7374.76	4816.85	379.79
					1150	760	40	4	3	3	5134.15	3360.49	249.02
					1380	912	48	60	59	59	90,812.20	59,527.32	4215.37
		Results according to [23]	[7:10, 7:20, 7:30, 7:40, 7:50, 8:00]	0.127	230	152	8	10	9	9	3800	2420	330
					460	304	16	10	9	9	6100	3940	410
					690	456	24	10	9	9	8400	5460	490
					920	608	32	10	9	9	10,700	6980	570
					1150	760	40	10	9	9	13,000	8500	650
					1380	912	48	60	59	59	91,800	60,120	4380
$r_{b\stackrel{}{\to }a}$	q = 4	Proposed method	[7:15, 7:28, 7:38, 7:46]	0.169	190	178	21	13	12	12	3680.26	3291.52	552.29
					380	355	41	10	9	9	4730.97	4301.94	624.84
					570	533	62	8	7	7	5304.77	4865.55	667.87
					760	711	83	90	89	89	76,778.71	70,757.42	9403.55
		Results according to [23]	[7:15, 7:30, 7:45, 8:00]	0.192	190	178	21	15	14	14	4800	4245	765
					380	355	41	15	14	14	7650	6900	1065
					570	533	62	15	14	14	10,500	9570	1380
					760	711	83	90	89	89	80,100	73,440	10,170
	q = 5	Proposed method	[7:12, 7:24, 7:35, 7:43, 7:49]	0.150	190	178	21	12	11	11	3562.43	3171.81	547.95
					380	355	41	11	10	10	5355.56	4854.49	722.28
					570	533	62	8	7	7	5414.95	4954.54	693.30
					760	711	83	6	5	5	5201.22	4783.91	645.97
					950	888	103	72	71	71	76,094.59	70,150.86	9191.68
		Results according to [23]	[7:12, 7:24, 7:36, 7:48, 8:00]	0.168	190	178	21	12	11	11	3840	3396	612
					380	355	41	12	11	11	6120	5520	852
					570	533	62	12	11	11	8400	7656	1104
					760	711	83	12	11	11	10,680	9792	1356
					950	888	103	72	71	71	77,760	71,496	9576

presents the combined results of the proposed method using MOPSO algorithm and the results from the literature. The results are displayed using the following parameters: the number of passengers left at a station, waiting time and amount of PWT. Each of these are displayed as a matrix with the number of columns representing the number of stops and the number of rows representing the number of services. Based on the input data and results in [23], the departure times from the terminals a and b, needed for comparison, are presented in Table 4, in three sets for the number of departures for route $r_{a\stackrel{}{\to }b}$ and two sets for the number of departures for route $r_{b\stackrel{}{\to }a}$.

With respect to the MOPSO algorithm, all experiments were performed using 200 particles (population size) and 500 generations. The detailed results are presented in Table 4.

In order to compare the results, the waiting time and amount of PWT for each solution are presented. As shown, the waiting time for the next service is shorter with the MOPSO algorithm although it is an uneven timetable. Hence, the amount of PWT based on (11) is better when using MOPSO. Although VOR is maximized (it is not displayed in Table 4, but it can be deduced by looking at the number of passengers left at the station and taking into consideration the fixed desired occupancy), the amount of PWT is not acceptable because of the high number of passengers left at the station and the long waiting time for consecutive service. From the obtained results, it can be concluded that passengers will choose another type of PT. In order to provide an example with a more user-oriented PT service, Experiment 1 is expanded as presented in Experiment 2.

5.2. Experiment 2

The second experiment is a bit more complex. All the parameters and conditions are the same as in Experiment 1 with two changes: the PT line is a train or tram line and the vehicles are passenger cars or coaches each having a capacity of 70. It is also assumed that a maximum of seven passenger cars can be connected. The aim of this experiment is to reduce the amount of PWT and maximize VOR at the same time. MOPSO was used to determine the optimal parameters. With respect to the MOPSO algorithm, all experiments were performed using 200 particles (population size) and 500 generations. The proposed optimal solution was found by the TOPSIS method. The obtained Pareto-optimal set and the proposed optimal solution is displayed in Figure 3, The detailed results are presented in

Table 5. Detailed results obtained in Experiment 2 using the proposed method.

				$1-{\mathit{\eta }}_{\mathit{av},\mathit{norm}}$	Vehicle Occupancy			Number of Free Seats			Number of Passengers Left at the Station			${\mathit{Z}}_{\mathit{av},\mathit{norm}}$	Waiting Time			Amount of PWT (Equation (11))
$r_{a\stackrel{}{\to }b}$	q = 4	Dep. time	[7:15, 7:25, 7:37, 7:46]	0.080	1	1	1	0	0	0	230	152	8	0.498	10	9	9	3348.39	2149.03	254.73
					1	1	0.52	0	0	236	40	136	0		12	11	0	1738.06	2386.84	209.68
		$d_o$	[70, 490, 70, 490]		1	1	1	0	0	0	270	288	8		9	8	8	3373.55	3158.13	229.26
					1	1	0.52	0	0	236	80	272	0		90	89	0	16635.48	30141.29	1572.58
	q = 5	Dep. time	[7:12, 7:21, 7:32, 7:40, 7:49]	0.064	1	1	1	0	0	0	230	152	8	0.057	9	8	8	3125.07	2001.04	247.84
					1	1	0.52	0	0	236	40	136	0		11	10	0	1729.53	2269.72	214.92
		$d_o$	[70, 490, 70, 490, 350]		1	1	1	0	0	0	270	288	8		8	7	7	3097.84	2866.7	220.31
					1	1	0.52	0	0	236	80	272	0		9	8	0	1775.07	3081.04	175.84
					1	1	1	0	0	0	310	424	8		72	71	71	30760.54	35592.32	1982.76
	q = 6	Dep. time	[7:10, 7:19, 7:29, 7:36, 7:45, 7:51]	0.080	1	1	1	0	0	0	230	152	8	0.040	9	8	8	3245.49	2073.29	267.91
					1	1	0.52	0	0	236	40	136	0		10	9	0	1669.51	2121.71	211.59
					1	1	1	0	0	0	270	288	8		7	6	6	2778.66	2549.20	204.11
		$d_o$	[70, 490, 70, 490, 70, 490]		1	1	0.52	0	0	236	80	272	0		9	8	0	1862.56	3133.54	190.43
					1	1	1	0	0	0	310	424	8		6	5	5	2621.71	3001.02	174.95
					1	1	0.52	0	0	236	120	408	0		60	59	0	14817.07	29050.24	1269.51
$r_{b\stackrel{}{\to }a}$	q = 4	Dep. time	[7:15, 7:25, 7:37, 7:46]	0.066	1	1	1	0	0	0	190	178	21	0.060	10	9	9	2808.60	2513.87	419.68
					1	1	0.63	0	0	154	30	194	0		12	11	0	1450.32	3208.65	251.61
		$d_o$	[70, 420, 70, 490]		1	1	1	0	0	0	220	371	20		9	8	8	2797.74	3999.48	368.71
					0.98	1	0.6	10	0	196	0	350	0		0	89	0	8177.42	38104.84	1887.10
	q = 5	Dep. time	[7:12, 7:22, 7:31, 7:41, 7:49]	0.055	1	1	1	0	0	0	190	178	21	0.061	10	9	9	2910.14	2595.88	443.11
					1	1	1	0	0	0	310	323	3		9	8	8	3699.12	3641.29	236.80
		$d_o$	[70, 140, 490, 70, 490]		1	1	0.57	0	0	212	80	306	0		10	9	0	1810.14	3875.88	233.11
					1	1	1	0	0	0	270	484	21		8	7	7	2968.11	4524.70	354.49
					1	1	0.60	0	0	194	40	469	0		72	71	0	10152.97	39642.32	1678.38

.

Comparing the results in Experiment 2 with those in Experiment 1, the waiting time for an uneven headway is acceptable when more than one passenger car or coach is used per service. The number of passengers left at the station is reduced, and it is assumed that passengers will wait for the next service.

6. Analysis and Discussion

In Experiment 1 (Table 4), the results indicate that by optimizing the timetable using the proposed method, an uneven timetable is obtained compared to the solution obtained by [23]. However, the amount of PWT is lower. In Experiment 2, it is assumed that the PT line is a train or tram line and that a maximum of seven coaches, each of capacity 70, can be used. The results obtained using the proposed method show a drastic improvement of the amount of PWT (i.e., a decrease in the value of $Z_{av,norm}$).

Table 6. Results obtained when vehicle occupancies of Experiment 2 are combined with corresponding departure times of Experiment 1.

	${\mathit{r}}_{\mathit{a}\stackrel{}{\to }\mathit{b}}$					${\mathit{r}}_{\mathit{b}\stackrel{}{\to }\mathit{a}}$
		Dep.Time	${\mathit{d}}_{\mathit{o}}$	$1-{\mathit{\eta }}_{\mathit{av},\mathit{norm}}$	${\mathit{Z}}_{\mathit{av},\mathit{norm}}$		Dep.Time	${\mathit{d}}_{\mathit{o}}$	$1-{\mathit{\eta }}_{\mathit{av},\mathit{norm}}$	${\mathit{Z}}_{\mathit{av},\mathit{norm}}$
q = 4	Proposed method	[7:15, 7:28, 7:38, 7:46]	70	0.080	0.050	Proposed method	[7:15, 7:28, 7:38, 7:46]	70	0.066	0.061
			490					420
			70					70
			490					490
	Results according to [23]	[7:15, 7:30, 7:45, 8:00]	70	0.080	0.064	Results according to [23]	[7:15, 7:30, 7:45, 8:00]	70	0.066	0.078
			490					420
			70					70
			490					490
q = 5	Proposed method	[7:12, 7:24, 7:35, 7:43, 7:49]	70	0.095	0.040	Proposed method	[7:12, 7:24, 7:35, 7:43, 7:49]	70	0.055	0.061
			490					140
			70					490
			490					70
			70					490
	Results according to [23]	[7:12, 7:24, 7:36, 7:48, 8:00]	70	0.095	0.048	Results according to [23]	[7:12, 7:24, 7:36, 7:48, 8:00]	70	0.055	0.073
			490					140
			70					490
			490					70
			70					490
q = 6	Proposed method	[7:10, 7:20, 7:30, 7:40, 7:47, 7:51]	70	0.080	0.040
			490
			70
			490
			70
			490
	Results according to [23]	[7:10, 7:20, 7:30, 7:40, 7:50, 8:00]	70	0.080	0.046
			490
			70
			490
			70
			490

shows the results that are obtained when the desired vehicle occupancy obtained in Experiment 2 is combined with the departure time obtained in Experiment 1. As was expected, both $(1-{\eta }_{av,norm})$ and $Z_{av,norm}$ values are decreased (implying an increase in VOR and decrease in the amount of PWT) compared to the solutions obtained by [23] in Experiment 1.

For example, for route $r_{a\stackrel{}{\to }b}$, when four departures and the possibility of using more passenger cars or coaches are considered, the value indicating the amount of PWT (0.050) is lower than that obtained by [23] $\left(0.064\right)$, while the value representing the vehicles’ occupancy ratio (0.080) is the same for both algorithms (Table 6). If the desired vehicle occupancy is fixed, i.e., $d_o=70$, the amount of PWT when using the proposed method is $0.142$, while it is $0.161$ when using the input given in [23] (Table 4). From these results, it can be concluded that it is more appropriate to use more passenger cars or coaches and an uneven timetable for the input data during the particular time period.

A qualitative analysis of the proposed timetable, with respect to the headway, is performed using the following indicators—the average headway (AH) and the expected waiting time (EWT), as recommended in [12]. AH and EWT of randomly arriving passengers, for all sets of departures and routes, are presented in

Table 7. Average headway and expected waiting time.

			Average Headway	Expected Waiting Time
$r_{a\stackrel{}{\to }b}$	q = 4	Proposed method – exp 1 / exp 2	10.33	5.37/5.24
		[23] – exp 1	15	7.5
	q = 5	Proposed method – exp. 1 / exp. 2	9.25	4.93/4.68
		[23] – exp 1	12	6
	q = 6	Proposed method – exp 1 / exp 2	8.2	4.45/4.23
		[23] – exp 1	10	5
$r_{b\stackrel{}{\to }a}$	q = 4	Proposed method – exp 1 / exp 2	10.33	5.37/5.24
		[23] – exp 1	15	7.5
	q = 5	Proposed method – exp 1 / exp 2	9.25	4.93/4.66
		[23] – exp 1	12	6

. For example, for route $r_{a\stackrel{}{\to }b}$ and six departures during the given time period, the expected PWTs, when using the proposed method, are 4.45 and 4.23 for Experiments 1 and 2 respectively, while the expected PWT is 5 when using the procedure in [23]. The presented qualitative analysis confirms that the timetable obtained using the proposed model and optimized using MOPSO is more appropriate compared to the timetable obtained in literature [23].

7. Conclusions

This paper presents a new model for determination of PT timetable. The model is formulated using a multiobjective optimization model to optimize VOR and PWT. Thus, this model takes into account the satisfaction of transport companies and passengers. The MOPSO algorithm is used in optimizing the model. The best solution in the Pareto-optimal set was found by the TOPSIS method. Practical implementation of the proposed model is presented using two numerical examples. PWT and VOR indices are used to present the performances of the proposed model in comparison to the similar results in the existing literature. Experiment 1 uses a simple passenger transportation line involving three sets of a number of departures (q = 4, q = 5 and q = 6) for route $r_{a\stackrel{}{\to }b}$ with four stations and two sets of number departures (q = 4 and q = 5) for route $r_{b\stackrel{}{\to }a}$. Experiment 2 has the same parameters and conditions as Experiment 1, and the additional assumptions that the PT line is a train or tram line and the vehicles are passenger cars or coaches. In both experiments, the proposed model using MOPSO algorithm shows better performances; i.e., shorter PWT and greater VOR. The case study based on the operation data from the existing literature shows that the proposed approach can reduce the average PWT by 10.54% for all sets with differing numbers of departures for two routes during the given time period. Based on the presented results, it can be concluded that the presented model, which uses the MOPSO algorithm to determine the optimal timetable, has an advantage in comparison to the existing models in scientific literature, which makes it suitable for scientists and practitioners in the field of PT. Further research will initially involve verification of the model using data from a real network scenario. We plan to modify the proposed model by including an extra constraint that ensures that the operation duration remains unchanged, especially if the number of vehicles or services is kept constant. We also plan to take into consideration boarding and alighting times, transfer stations and scheduling several lines.