Probabilistic field approach for motorway driving risk assessment

We present an approach to assess the risk taken by on-road vehicles within the framework of artificial field theory, envisioned for safety analysis and design of driving support/automation applications. Here, any obstacle (neighboring entity on the road) to the subject vehicle is treated as a finite scalar risk field that is formulated in the predicted configuration space of the subject vehicle. The driving risk estimate is the strength of the risk field at the subject vehicle’s future location. This risk field is formulated as the product of two factors: collision probability and expected crash energy. The collision probability with neighboring vehicles is estimated based on probabilistic motion predictions. The risk can be assessed for a single time step or over multiple future time steps, depending on the required temporal resolution of the estimates. We verified the single step approach in three near-crash situations from a naturalistic dataset and in cut-in and hard-braking scenarios with simulation and showed the application of the multi-step approach in selecting the safest path in a lane-drop section. The risk descriptions from the proposed approach qualitatively reflect the narration of the situation and are in general consistent with Time To Collision. Compared to current surrogate measures of safety, the proposed risk estimate provides a better basis to assess the driving safety of an individual vehicle by considering the uncertainty over the future ambient traffic state and magnitude of expected crash consequences. The proposed driving risk model can be used as a component of intelligent vehicle safety applications and as a comprehensive surrogate measure for assessing traffic safety.


Introduction
Traffic safety has attracted increasing research attention, in particular in the transition from human-driven vehicles to automated vehicles. In general, safety can be defined as the absence of unreasonable risk. Hence, driving risk models are key components of various integrated vehicle driving applications such as collision warning/avoidance systems. Modelling dynamic driving risk entails a detailed description of the subject vehicle and its environment. Fortunately, safety analysts now have the means to scrutinise driving at a far more detailed level, thanks to the high resolution driving data provided by modern sensing and communication technologies. In this context, we explore the possibility of developing an assessment method for driving safety (converse risk). Such a method can be used to assess the risk of human driving, and to evaluate (and design) proactive safety systems and advanced vehicle control systems. Towards this, we survey the relevant literature from the safety analysis domain and propose a method for modelling driving risk.

Objective and contribution
Using the paradigm of artificial field theory, this paper puts forward an approach to assess the driving risk for on-road vehicles named Probabilistic Driving Risk Field (PDRF). The approach takes into account the probability of motion predictions of other road users. For an encounter between the subject vehicle and a road obstacle, driving risk is defined as the consequence to which the subject exposes itself by maintaining its present kinematic state in spite of the uncertainty that the obstacle's future motion may (or may not) lead them to collide. The driving risk estimate constitutes a crash severity term and a collision probability term (ISO 26262-1, 2018;Aven et al., 2018). The subject and neighboring vehicle's possible positions and associated probabilities at discrete future time steps are predicted. The collision probability is modelled as the likelihood of overlap in the predicted spatial configurations of the subject and neighbor vehicle, under the assumption of motion prediction. The risk can be estimated for a single time step or over multiple future time steps, depending on the required temporal resolution. We validate the single-step approach in near-crash scenarios from naturalistic driving data and simulations of typical risky scenarios. Additionally, we apply multi-step risk to evaluate the risk along four trajectory options to the subject vehicle while approaching a lane drop. This serves to evaluate risk for alternative plans generated by an advanced vehicle controller.
In the remainder, Section 2 formulates the risk model in the single-step approach while Section 3 extends it to multiple steps. Section 4 presents the simulation experiments and empirical validation of the proposed model, while Section 5 reflects on the advantages and limitations of the model. Section 6 concludes the paper.

Modelling driving risk in the single-step approach
We first present the models of the risk fields in the single-step approach. Any obstacle that the subject vehicle encounters while driving on the road is formulated as a risk field in the configuration space of the subject vehicle. The strength of the risk field at a spatial point around the obstacle is formulated as the product of the crash severity and collision probability with the subject vehicle.
The road space is considered as a flat Euclidian plane consisting of road surface markings. In this work, we model the driving risk along a straight road stretch, and hence lane markings are treated as straight lines. The axis X is in the direction of the vehicle movement and is aligned along the outer lane marking, and the axis Y is in the direction perpendicular (counter-clockwise) to the  Fig. 1a. We define two types of on-road obstacles: road boundary objects existing outside the driving lane such as medians, roadside barriers; and movable obstacles such as vehicles (Wang et al., 2016). All the road boundary objects are treated as an immovable entity, and medians and barriers are modelled as finite line segments parallel to the X axis Vehicles are modelled as a rectangle with its length denoted as L and width denoted as W, and possess a finite physical mass.
In our Probabilistic Driving Risk Field (PDRF) framework, the road boundary objects are modelled as potential PDRF and vehicles are modelled as kinetic PDRF. The risk estimate is the strength of a PDRF at the location of the subject vehicle, s. Accordingly, at a given moment, s takes a total risk R which is comprised of potential risk (total potential PDRF strength) R B from multiple road boundary objects and kinetic risk (total kinetic PDRF strength) R N from multiple dynamic road objects as follows: where b denotes an individual road boundary object and n denotes the individual vehicle. Based on the superposition property of fields, the total risk model of Eq. (1) combines risk posed by multiple road obstacles into a single measure. This implies that the risk posed by an obstacle is assumed to be independent of the risk posed by another obstacle. The risk is formulated as the product of crash severity and collision probability. Vehicle crashes are not perfectly elastic (collision with no loss of net kinetic energy); and some part of the kinetic energy is dissipated as thermal energy, sound energy, and material deformation. Given a crash, the portion of dissipated energy that is spent on deforming s is termed as expected crash energy. We use expected crash energy as the approximation of crash severity. The collision probability will be specifically formulated for fixed and moving obstacle types.

Road boundary object as a risk field
The potential risk taken by s due to a fixed road boundary object b is modelled as follows: where R b s , denotes the strength of the PDRF due to the road boundary object b M ; denotes the mass of s r ; s b , and r L are vectors that denote the shortest distance between b and s (the centre of gravity) and lane center respectively; V s b , denotes the velocity of s along r s b , ; The potential risk R b s , can be separated into two terms: a crash severity term and a collision probability term. The crash severity term, 2 denotes the expected crash energy scaled by the parameter k, with range [0-1], representing the rigidity of the road boundary object. Accordingly, = k 1 when b is immovable with infinite mass; and = k 0 when b is very compliant and effectively dissipates the crash energy (Zou et al., 2014). Among road boundary objects, the cable barrier, in general, dissipates the highest amount of energy; followed by the guardrail, the concrete median and any fixed objects such as trees and roadside poles.
The collision probability term e r D , . This is intuitive; a road object further away offers more possibility for the driver to avoid the collision (Zou et al., 2014). . The influence of R b is restricted to a finite area by two means. First, we set . Second, we set = D |r |/7 L meaning that collision probability term attains a marginal value (0.00091 0.001) at the lane centre. To illustrate the sensitivity of the risk, let us consider the situation illustrated in Fig. 1a where the lane boundary b aligns with the X axis and is positioned at an offset from the lane boundary, denotes the angle (counterclockwise) between V s and X axis. Fig. 1b illustrates the sensitivity of the risk to the vehicle heading angle. It can be seen that risk is highest when the vehicles is moving straight onto the barrier ( =°90 ). Fig. 1c illustrates the sensitivity of the potential risk towards the parameter k and the offset of b.

Modelling the neighboring vehicle as a risk field
The kinetic risk taken by s due to a neighboring vehicle n is modelled as follows 2 describes the portion of expected crash energy that has to be absorbed by s if it collides with n. The collision process is assumed to be inelastic where both the vehicles move together after their first contact. Here, the division of crash energy is inversely proportional to the individual mass. Therefore a lighter vehicle will dissipate more energy than a heavier vehicle. The second term p n s ( , ) with range [0, 1] denotes the crash probability. According to the model of Eq. (3), when = p n s R ( , ) 1, n s , attains a finite maximum = R M V 0.5 | | n s s n , , 2 , i.e. the expected crash energy. A vehicle, unlike a road boundary object, is a movable entity and therefore the collision probability term in Eq. (4) takes account of the dynamic nature of the interaction. We operationalise the risk definition as follows: risk is the consequence of s maintaining its state (single motion prediction assuming an acceleration signal, A = 0 m/s 2 ), despite the unknown motion of n (a continuous range of acceleration signals A A A { , } min max ). Earlier studies have used constant acceleration heuristics to calculate TTC along a single motion prediction ; but here we consider a range of motions possible to n.

Vehicles dynamics model for motion prediction
The dynamic state of a vehicle i is described in global coordinates by the position of its centre of mass = P X Y [ , ] i i Xi Yi T with V Xi and V Yi respectively denoting the longitudinal and lateral components of the velocity. Vehicle yaw rotation is not modelled explicitly. Instead the vehicle motion direction (heading) follows from the vector V. A vehicle is treated as a dynamic system whose state is manipulated by its input which is acceleration Xi Yi T , with A Xi and A Yi respectively denoting the longitudinal and lateral components of the acceleration. In order to predict the evolution of the vehicle state, the dynamic behaviour of the vehicle is modelled as the following differential equation: The vehicle motion is subject to physical constraints. We emulate them in the vehicle model as a set of constraints: X representing the non-holonomic behaviour of the vehicle. This condition assumes that the vehicle heading angle = ( ) is bounded as | | 10 , during motorway driving; representing the strictly forward movement; • A A A X X X min max representing the feasible acceleration range that is restricted by the engine power and brake torque limitations.

Collision probability with the neighbor
Given the above future motion of s and n, we now estimate the probability of collision. In this work, we estimate the collision probability at a single future time instant; therefore the collision probability is merely related to the chances of spatial overlap. Accordingly, we use the prediction time step as an estimation parameter, and p n s ( , | )denotes the probability of collision (spatial overlap) between n and s at a future time + t 0 . To estimate p n s ( , | ), we only use the predicted position of s and the range of predicted positions of n at + t 0 . This predicted approach can be viewed as the simplest case in the tree of possible paths: with a single time step which is also the prediction horizon. Note that this collision probability estimate is conceptually different from SMoS such as TTC, which approximates the probability of a spatial and temporal overlap. A human driver requires a finite duration commonly referred to as the reaction time (Treiber et al., 2006) to sense, perceive, decide and act to stimuli. Therefore, s is unaware and unresponsive to the action of n within this reaction time. It must be noted that s cannot perceive kinetic risk estimated for shorter than its reaction time.
The probability of motion predictions is attributed to the underlying acceleration signal. The probability functions of acceleration variability can be estimated by treating acceleration signals as a random variable (Wagner et al., 2015). We assume acceleration variability to follow a normal distribution (Ko et al., 2010). Accordingly, the expected motion (i.e. the vehicle motion, assuming the acceleration signal equal to the mean of the acceleration variability distribution, typically 0 m/s 2 ) is associated with the highest probability. The parameters of this distribution can be estimated as follows: where µ A denotes the mean acceleration; A denotes the standard deviation of the acceleration (acceleration noise) and T denotes the sampling time duration. The acceleration noise when a vehicle is stopped in traffic is zero, which might distort its estimated value; therefore the acceleration noise should be measured only while the vehicle is moving (Jones and Potts, 1962).
The absolute probability for a continuous random variable to take a particular value is 0, and it can only be determined for a particular range, which will be defined in the next section by accounting for the vehicle geometry. To illustrate the concept and evaluate the mathematical properties, we express each vehicle as an infinitesimally small point mass. Here, the value of the acceleration variability distribution for a particular value of acceleration is interpreted as the relative likelihood of occurrence. We tentatively consider "collision likelihood" defined as the relative likelihood of n applying acceleration: . Considering A X n , and A Y n , as random variables, the collision likelihood p L can be defined as follows: where is the probability density function, and its parameters µ denotes the mean and denotes the standard deviation of the distribution; = X X X s n and = Y Y Y s n denotes the relative spacing, and , denotes the relative velocity in longitudinal and lateral directions. By substituting for p n s ( , ) in Eq. (4) with p L in Eq. (6), we obtain a specific form of R n s , . The contours in Fig. 2 represent the relative positions of s around n where the collision likelihood is the same. As described in the introduction, the acceleration noise can be impacted by road, traffic state and driver-related factors. It can be seen that the PDRF approach relates an increase in acceleration noise to an increased region of risk; the region of high collision likelihood widens with an increase in the acceleration noise (see Fig. 2(b) and (c)). Furthermore, it can be seen that the collision likelihood is sensitive to the parameter (see Fig. 2(d)).

Incorporating the vehicle geometry and motion constraints in the collision probability
We introduced the kinetic PDRF field based on the estimates of collision likelihood between point mass vehicles, treating the acceleration signal as a continuous random variable. However, the vehicle possesses a geometry and a crash is the overlap in geometric boundaries of vehicles. In this section, we specify the collision probability of the vehicle model specified in Section 2.2.1 (rectangular geometry and constrained motion). Imposing the motion constraints in Section 2.2.1, the boundary of the reachable state of n at time + t 0 can be represented as quadrilateral polygon Q. Using the predicted position of s at + t 0 along with the geometry of the two vehicles, we define the potential collision zone: Z. Thereby, the collision probability is non zero if there exists an overlap in Q and Z. The overlapping region is another polygon denoted by O defined in the spatial domain. The calculation of Q Z , and O are provided as an implementable algorithm in the appendix A. The overlap O in the spatial domain has to be converted to acceleration domain denoted A O for probability estimation by using the following relation.
where X Y , c c denotes the corner positions of O, and A A , X c Y c denotes the corresponding corners of A O . Eq. 3.7 specifies a linear relationship between acceleration coordinates and the spatial coordinates, and therefore A O is also a quadrilateral polygon. Then the collision probability can be obtained by integrating the joint acceleration variability distribution over A O as follows:

Estimating risk over multiple time steps
In the previous section, we described the procedure to estimate risk based on collision probability for one fixed time step. Such a risk estimate is intended to assess the driving safety of human drivers from observed vehicle interactions; and to be used as a trigger for collision warning/avoidance systems. However, advanced vehicle control systems plan the vehicle trajectory over longer time horizons. Therefore, risk estimates over a sequence of future time steps are necessary to create and evaluate such trajectory plans. The single-step PDRF approach is inadequate for this purpose, and therefore, we extend the definition of the PDRF model. The total risk taken by s at each time step t k is R t ( ) T k which is defined as follows:  where R t ( ) B T k denotes the potential risk (total potential PDRF strength) at time t k due to the presence of road boundaries; and R t ( ) N T k denotes the kinetic risk (total kinetic PDRF strength) at time t k due to the presence of neighboring moving objects. R t ( ) b T k is the potential risk due to the individual road boundary object b and is equivalent to R b (modelled in (3.3)) calculated at t k using the expected conditions at the time t k . For path planning purposes, in addition to the road boundary objects that can inflict physical crash energy, it might be of interest to formulate virtual road objects such as lane marking as artificial risk fields. This can be achieved by adjusting the value of the scaling parameter k to match the subjective level of risk attached to the lane marking type.
In order to extend the model (2) towards generating risk estimates over multiple time steps as in Eq. (9), three questions must be answered: (1) How to predict the possible motions of over multiple time steps? (2) How to determine the probability of these motion predictions? (3) How to estimate the collision probability under a probabilistic setting? The answer to the third question can be obtained from the work of Saunier et al. (2010). They illustrated a procedure to identify the collision points along a pair of motion predictions, and thereby to express collision probability as the discrete sum over a finite number of collision points. Questions 1 and 2 can only be answered by considering the envisaged application.
The multi-step prediction is envisaged to estimate the safety of trajectory plans. Trajectory planning is done when the subject vehicle s is operated by an advanced vehicle control system. The motion possibilities of neighboring vehicles are predicted by a scheme which employs two functions: (1) an acceleration plan which is defined as the discrete time series of expected acceleration signals that span a finite prediction horizon; (2) the spread of the acceleration signal at each time step along the discrete time series. The acceleration signal at each time step is considered as a random variable whose variability distribution is parameterised by the expected acceleration and its spread. The horizon of an acceleration plan determines the length of the trajectory prediction and spread determines the spatial spread of the predicted trajectory set. In this work, the multi-step approach is used to assess the risk of trajectory plans. However, an advanced vehicle control system can use this approach to evaluate the safety of its candidate trajectory plans. For this, the system should predict the acceleration plan of n (using manoeuvre prediction algorithms) and the spread in acceleration signals (using on-board measuring and estimation systems).
We answer the first question by describing the approach to predict the motion of n over multiple time steps. According to the risk definition, if s maintains its trajectory plan and would be exposed to the consequence of the uncertain motion of n (a finite set of trajectory predictions).The trajectory of a road user i is denoted as is the single trajectory plan of s in time discrete form, with k as the prediction time instance and H as the total number of prediction time steps; denotes the finite set of possible trajectories of n, where m denotes the trajectory instance and M denotes the total number of predictions. The motion possibilities of n are predicted as a tree of possible paths (Jansson, 2005). The motion of the neighboring vehicle is modelled as a wiener process evolving at discrete time steps governed by the acceleration input which is treated as a random variable. To predict one trajectory instance, a pair of the acceleration signals A A ( , ) X Y is assigned to n at the beginning of each time step, and thereafter it is propagated using the assigned accelerations during the time step. In this study, we use a set of 25 unique pairs of acceleration inputs, i.e. all possible combinations among the 5 values (with an interval of 1 m/s 2 ) of A X and A Y . At each time step, all the previous end states are assigned this set of accelerations and they are further updated. This process is done at each time step until the end of the prediction horizon (H = 4 s in this work) and results in a set of 390625 trajectories. After that, we impose motion constraints to the predicted trajectory set to eliminate the infeasible trajectories. Now we answer the second question by describing the method to assign a probability to each trajectory instance within the finite set of predictions using the acceleration probability distribution. Since the trajectory predictions are finite, we convert the continuous probability density function to a discrete probability function. The discretisation interval is 1 m/s 2 , corresponding to the acceleration interval used for motion prediction. Thereby, the probability assigned to each trajectory prediction is the scalar product of individual probabilities of the underlying sequence of acceleration signals.
. The acceleration variability distributions follow a discrete normal distribution as follows: . We now deploy the approach proposed by Saunier and Sayed (2009) to estimate the crash probability for a finite set of trajectory predictions. Towards this, we first define collision event as the function Proximity A B ( , ) defined for a given spacing d and spacing threshold as follows: In order to calculate the collision probability, we need to identify those trajectories that will lead to a collision and sum up their probability. Upon checking for the condition = T T Proximity ( n m H , is the first time instant satisfying the proximity condition. Let g k be a function defined over CP specific to t k that returns all the trajectory predictions of n those lead to the collision. The crash probability can then be estimated as follows: where C k is the total number of collision points at t k . Now R t ( ) n T k is a time series of R n calculated at t k using the expected conditions at time t k and p s n t ( , | ) k

Example model applications
In this section, we illustrate the applicability of the risk assessment approach and evaluate its performance. First, we apply the single step PDRF model to generate risk estimates along vehicle trajectories observed in three near-crash events selected from a naturalistic driving dataset (Neale et al., 2005). This will reveal if the risk descriptions match the event narration and provide insights into the safety-critical moments during the event. As seen in Section 2 (see Fig. 2(d)) the risk measure is influenced by the prediction time step . Therefore, we will use this experiment to evaluate the sensitivity of . Secondly, we evaluate the single-step PDRF, by simulating several variations of two typical hazardous scenarios: cut-in and hard braking. This test will reveal the accuracy and error rate of PDRF metric, in these two scenarios. In the next step, we will apply the multi-step approach to estimate the risk of four possible trajectories that the subject could pursue while approaching a typical motorway lane drop section. This will reveal if the estimates are reasonable and allow comparison to a standard safety indicator: generalised TTC. Besides, we will use this experiment to illustrate the application of the motion prediction scheme to calculate the generalised TTC proposed by Saunier et al. (2010).

Risk assessment of near-crash events
The sample trajectories are obtained from the public dataset provided by the 100 vehicles naturalistic study (Neale et al., 2005). The study employed volunteering drivers with instrumented vehicles to collect large-scale naturalistic driving data. Vehicles were instrumented with data collection equipment including video cameras, front and rear radars; GPS positioning system and vehicle motion sensors. This database provides a detailed record of 761 near-crash situations in terms of annotated event video and a textual narration of the incident. The dataset includes situations labelled as "near-crash", meaning "a conflict situation requiring a rapid, severe evasive manoeuvre to avoid a crash" (Neale et al., 2005). The presence of an evasive action indicates that the driver considered the situation to be risky, and therefore these situations can be regarded as examples of conflicts characterising "true risk". The dataset provides time-stamped vehicle states of the subject vehicle and its spacing (range) and relative longitudinal velocity (range rate) with the front and rear vehicles (maximum 7 vehicles). However, the dataset does not contain the lateral kinematics of neighboring vehicle, and hence its trajectory cannot be reproduced in detail. We selected three near-crash situations with solely the ego vehicle, or those in which it was clear from the textual narration that the neighboring vehicle did not perform any lateral maneuver. Fig. 4 illustrates the three traffic scenarios: (1) a road boundary crash that was evaded by corrective steering; (2) a rear end near-crash that was evaded by swerving to the left; (3) a rear-end crash that was evaded by braking. The trajectory samples were extracted via the following steps: (1) identify the safety critical neighboring vehicle in the encounter, based on the textual narration of the event; (2) extract the trajectory of the subject and the safety critical neighboring vehicle; (3) estimate the lateral velocity (longitudinal acceleration) from the difference in lateral position (longitudinal velocity) between two consecutive time records. The lateral velocity Fig. 4. Illustration of the near-crash scenarios. F.A. Mullakkal-Babu, et al. Transportation Research Part C 118 (2020) 102716 and longitudinal accelerations were found to be highly noisy, and therefore we applied a moving average filter with a time span of 1 s.

Road boundary crash avoided by corrective steering
We consider the near-crash event during trip number 8299 in the study dataset. The subject vehicle is in the entrance/exit only lane, trying to get into the left lane. The subject inadvertently drifts towards the concrete barrier on the right side of the road. She applied corrective steering to bring the vehicle back to the lane centre. We calculate the total risk R taken by the subject vehicle according to Eq. (2), = R R b where b is the concrete barrier. Since the offset of the concrete barrier from the lane boundary is unknown, we consider the barrier to be at an offset of 0 m. We set the value of scaling factor k to 0.61 (the odds of injury ratio with a concrete barrier wall measured by Zou et al. (2014). Additionally, we calculate the Time to Lane Crossing (TLC, the ratio of lateral spacing to the lateral velocity of the subject in the direction of the barrier). The lateral dynamics of the subject are described in Fig. 5 (a), (f). Both TLC and PDRF risk descriptions qualitatively reflect the event narration in the dataset. During the initial phase of the encounter (time 0-2 s), both PDRF and TLC risk estimates remained marginal. Thereafter, the subject drifted towards the right barrier which increased the chances of crashing onto the right boundary; this unsafe development is indicated as a rise in the PDRF risk description (see Fig. 5 (d) at 2 s); and as a descent in the TLC description. Realising the impending hazard, the subject applied corrective steering (see Fig. 5 (d) at 3.2 s) to bring the vehicle back to the lane centre; the success of evasive action is indicated as a descent in the PDRF risk description from a maximum (see Fig. 5 (d) at 3.2 s) to a marginal value (see Fig. 5(d) at 3.6 s); and as an ascend in the TLC description. Fig. 5(c) describes the variation of the crash severity term and Fig. 5(c) describes the variation of the crash probability term in PDRF risk as given in Eq. (4). As expected from its formulation, the crash probability term correlates with the lateral position (see Fig. 5 (a) and Fig. 5 (c)). It can be seen that the PDRF approach provides continuous risk estimates with relatively smooth variations, and therefore it could clearly mark critical moments during the encounter. TLC goes to infinity as soon as the vehicle is no longer moving to the right (at 3.8 s), and therefore does not capture the risk due proximity of the rigid barrier. The risk was first seen to appear at 2 s (in Fig. 5 (d)). This provides an approximate time point to investigate driver actions or circumstances that lead to the event. Secondly, the risk estimate attained the maximum at 3.2 s (see Fig. 5 (d)). The maximum risk provides a representative risk measure for the entire encounter.

Rear crash avoided by swerving
We consider the near-crash event during trip number 8450 in the study database. The preceding vehicle brakes and slows down on the interstate. The subject vehicle is fairly far behind the preceding vehicle, but the preceding vehicle decelerates at a faster rate than the subject driver is initially aware of, so the subject has to brake and steer left to avoid hitting the preceding vehicle in the rear. We calculate the total risk R taken by the subject. According to Eq. (2), = R R n , where n is the preceding vehicle. The trajectory of n was not available long enough to estimate the distribution of acceleration. Therefore, we used a standard set of parameters to define the acceleration variability of n as A ( | 0, 0.7) . The evasive manoeuvre undertaken by the subject can be seen from its lateral movement in Fig. 6 (a); and the applied braking in Fig. 6 (f). Both TTC and PDRF risk estimate qualitatively reflected the event narration. During the initial phase of the encounter (time 0-6 s) the subject vehicle was far behind suggesting a safe following. During this phase, both PDRF risk (see Fig. 6 (d)) and TTC (see Fig. 6 (e)) estimate depicted a marginal risk. Thereafter, the preceding vehicle suddenly decelerates at a fast rate; this unsafe development is described as a steep descent in TTC (8-10 s) and as a steep and temporally adjacent rise in the PDRF risk descriptions ( = 4, 3 and 2 s).
It can be seen that the PDRF approach provided continuous risk estimates with relatively smooth variations (see Fig. 6 (d)) within a finite range; and the crash probability and crash severity varied independently (see Fig. 6 (b) and (c)). As expected, the prediction time step ( ), influenced the risk estimates, more specifically the collision probability term. Therefore the PDRF risk descriptions differed temporally and quantitatively (see Fig. 6 (d)).
However, evaluating the multiple descriptions of risk (see Fig. 6 (d)) and crash probability (see Fig. 6 (c)) provides information about the evolution of risk during the encounter. The peak in crash probability with = 2 s indicates high chances of an imminent crash at 2 s. Moreover, subsequent and comparable risk peaks ( = 2 and 3 s) around 9 s suggest that the braking, which started around 8 s is not sufficient to avoid the danger. In combination, these observations indicate a growing urgency for effective evasive  F.A. Mullakkal-Babu, et al. Transportation Research Part C 118 (2020) 102716 action. This could explain why the subject swerved to its left (at approximately 9 s) instead of merely braking. Moreover, the absence of crash probability with = 1 s indicates the collision was successfully avoided at least until the end of the observation. Notably, the PDRF risk with = 3 s attained the highest estimate and its moment of maximum risk is closer to the moment when the subject begins the evasive manoeuvre: swerve with hard braking around 9 s (see Fig. 6 (a), (f)).

Rear crash avoided by braking
We consider the near-crash event during trip number 8427 in the study database. The subject vehicle is following the lead vehicle fairly closely on a 2-lane road when the preceding vehicle slows down to stop. Subject vehicle brakes to avoid hitting the preceding vehicle in the rear. Similar to the previous experiment, we calculate R n and TTC estimates. Both TTC and PDRF risk estimates qualitatively reflect the event narration. During the initial phase (0-4 s), the subject was closely following the preceding vehicle indicating an unsafe interaction. During this phase, both PDRF risk (see Fig. 7 (d), = 4 s) and TTC (see Fig. 7 (e)) indicate the existence of risk. Thereafter the preceding vehicle slows down to stop; this unsafe development is described as a gradual descent in TTC (time 4-6 s) and as a temporally adjacent rise in the PDRF risk descriptions ( = 4, 3 and 2 s).
Evaluating the multiple descriptions of risk (see Fig. 7 (d)) and crash probability (see Fig. 7 (c)) reveals how this situation differs from the previous one in terms of risk evolution. The low peak in crash probability with = 2 s indicates marginal chances of an imminent crash (i.e. at the next 2 s). Moreover, a high-risk peak with = 3 s around 5.6 s and a subsequent lower risk peak with = 2 s suggest that the evasive braking, which started around 5 s, is sufficient to evade the danger. In combination, these observations indicate that the braking was effective to evade the danger. The absence of risk with = 1 s indicates the collision was successfully avoided at least until the end of the observation. Similar to the previous example, the PDRF risk with = 3 s attained the highest estimate and its moment of maximum risk (see Fig. 7 (c)) is closer to the moment when the subject begins the evasive braking around 5 s (see Fig. 7 (f)). It can be seen that the PDRF risk model could qualitatively reflect the event narration, and its risk description was consistent with that by TTC. The prediction time step influenced the PDRF risk descriptions. When evaluated together, the risk descriptions with distinct prediction time steps could provide information about the risk evolution of the encounter. However, each of these descriptions differ in their qualitative properties. For example, a peak in crash probability with = 1 s cannot be used to detect near-crashes; as it did not appear in the two near-crash examples. Similarly, a peak in crash probability with = 4 s cannot be  F.A. Mullakkal-Babu, et al. Transportation Research Part C 118 (2020) 102716 regarded as an accurate indicator of near-crash; such peaks occurred multiple times during a single encounter. Notably, in both the examples, the PDRF risk estimate with = 3 s yielded a single peak; attained the maximum value; and was temporally closest to the moment when the driver initiated the evasive manoeuvre.

Simulation based validation
To further validate the approach, we perform simulations of two typical safety-critical scenarios, and evaluate the accuracy of single-step PDRF ( = 3 s) in identifying the risky situations, in comparison to TTC. The simulated scenarios include both crash and non-crash events. TTC < 3 s and PDRF > 0 Joules are treated as positive risk flags, and crash between ego and neighbor is treated as TRUE risk. We quantify the accuracy of risk estimates in terms of the True Positive, True Negative, False Positive and False Negative risk descriptions.

Performance evaluation in cut-in events
In this scenario, ego vehicle travels on the left lane and neighbor vehicle travels on the right lane with a constant longitudinal velocity. As shown in Fig. 8a, the neighbor is 15 m ahead of the ego at t = 0 s, and cuts in onto the left lane at t = 6 s, with a constant lateral velocity of 1 m/s. This scenario is replicated by varying the initial velocity of ego between [5-30 m/s] and that of neighbor between [5-30 m/s] with an increment of 1 m/s resulting in 676 (26x26) scenario instances. Among the simulated scenarios, 49 instances result in a crash, as shown in Fig. 8a. It can be seen in Fig. 8b that PDRF correctly classified the scenario instances. In contrast, the false negative (see black region in 8c) reveals that TTC fails to capture the true risk in some scenarios, about 50% of the scenarios using a threshold of TTC<3 s. The acceleration distribution of the neighbor vehicle is set as

Performance evaluation in hard braking events
In this scenario, the ego vehicle travels on a lane and a leading vehicle travels in the front on the same lane starting with a constant longitudinal velocity, as shown in Fig. 9a. At t = 6 s, the lead vehicle abruptly applies a constant deceleration of 5 m/s 2 to stop. We vary 3 parameters of this scenario. The initial spacing between the two vehicles is set as (20 m, 40 m, 60 m, 80 m). The initial velocity of the preceding vehicle is varied between [5-30 m/s] with an interval of 1 m/s and the maximum velocity of the two vehicles corresponding to each of the initial spacing is 10 m/s, 16 m/s, 23 m/s and 30 m/s respectively. This results in 1217 simulations (36, 144, 361 and 676 simulations for different initial spacing respectively). To avoid crashes before hard braking, the minimum initial time headway between the two vehicles is set to 3 s. Therefore the initial velocity of ego is varied between a minimum of 5 m/s to a maximum of 30 m/s. Comparing Fig. 9c and 9d, it can be seen both PDRF and TTC result in few false-positive cases, particularly in those scenarios in which the relative velocity between the ego and preceding vehicle is relatively small. These  F.A. Mullakkal-Babu, et al. Transportation Research Part C 118 (2020) 102716 false positives represent high risk interactions leading to near-miss events. The acceleration distribution of the preceding vehicle is set as

Summary of simulation results
Compared to TTC, it can be seen that the number of false risk values is lower and true risk values are higher in PDRF risk estimate (see Table 1). This shows that PDRF detects the dangerous situations more efficiently and accurately than TTC, consistently in all the scenarios.
The simulation experiments demonstrate that PDRF correctly captures the crash and non-crash events in the cut-in scenario. In contrast, TTC resulted in 24 false-negative cases. These cases were sideswipe crashes when the cut-in vehicle approached the ego vehicle in the lateral direction.
In hard braking scenarios, both TTC and PDRF flagged positive risk in true crash events. In scenarios when the ego approached from an initial spacing > 40 m, the positive values of PDRF estimates were more reliable than TTC, as they contain lesser false positive values. However when the scenario with initial spacing of 20 m, both PDRF and TTC had the same error rate.
In all the true positive cases the value of TTC is 0. In contrast, the maximum PDRF values varied among the true positive cases depending on the expected crash severity. This shows the ability of the PDRF to scale the risk according to the severity.

Multi-step risk estimation of path plans
In this section, we apply the multi-step PDRF to estimate the risk of four trajectory plans while approaching a typical lane drop section. The situation involves four vehicles: s (the subject vehicle), n n n 1, 2, 3 (neighboring vehicles). The initial states are shown in Fig. 10. Here, the controller of s makes four trajectory plans: T1, T2, T3 (lane change trajectories) and T4 (forward trajectory). All the lane change trajectory plans span a duration 4 s and imply a constant lateral velocity of 1 m/s. Even though these trajectory plans are   F.A. Mullakkal-Babu, et al. Transportation Research Part C 118 (2020) 102716 similar in terms of lateral movement, they differ in the prescribed longitudinal dynamics: T1 represents a constant velocity of 20 m/s; T2 represents a decreasing velocity with a constant deceleration of −2 m/s 2 , T3 represents a decreasing velocity with a milder constant deceleration of −0.7 m/s 2 . T4 represents moving forward at 20 m/s in the original right lane. In this experiment, we set the trajectory plans A t ( ) X k and A t ( ) Y k of neighboring vehicles as shown in Fig. 11. The acceleration plan of n3 implies that it will continue to accelerate at 1 m/s 2 (see Fig. 11 (c)); n1 and n2 according to the Intelligent Driver Model (see Fig. 11 (a)); and n2 will move to the middle lane (see Fig. 11 (b)). The acceleration variability distribution is set as follows: According to the PDRF model, the risk of a crash between s and n1 is highest along T2 (lane change with hard braking); both in terms of collision probability and crash severity (see Fig. 12 (a)). Whereas, in the encounter between s and n2, the crash probability is highest along T1 (lane changing at constant velocity) (see Fig. 12 (b)). The risk description indicates a significantly high risk along  T2; and lowest along T4 (see Fig. 12 (d)). This means that objectively the safest trajectory is T4; however, if the subject decides to change the lane, it should avoid T2, and choose between T1 and T3. Fig. 12 (d) shows that even though the maximum risk values along T1 and T3 are comparable (T3 slightly lesser than T1), the danger is more imminent along T3 (at 3 s) than along T1 (at 4 s). The experiment shows the applicability of multi-step PDRF to generate ex-ante risk estimates in order to differentiate the trajectory plans. To compare our approach, we describe an additional metric generalised TTC proposed by Saunier et al. (2010). By this generalised definition of conflict, TTC is the time required for two vehicles to collide, following the predicted trajectories. The generalised TTC for the set of trajectory predictions is defined as follows:  Fig. 13 describes the generalised TTC values along the four trajectory plans. p TTC is the aggregated collision probability for a given t k , according to 3.2.5, which can also be interpreted as the TTC. Fig. 13 shows the generalised TTC, with H = 4, along the four trajectory plans. The trend shown by the generalised TTC values reflects the temporal proximity of a crash along the trajectories. Both generalised TTC and PDRF estimates suggest that T4 is the safest trajectory. Notably, the riskiest trajectory is T3 as per generalised TTC estimates (Fig. 13), whereas it is T2 as per PDRF (see Fig. 12 d). This is due to the absence of severity dimension in generalised TTC. It can be seen that the high crash severity is the factor that differentiates T2 for encounters with n1 and n2 (see Fig. 12 (a) and See Fig. 12 (b)). Even though generalised TTC is an effective and simple approach to detect critical interactions, multi-step PDRF risk contains the crash severity information, which is vital to differentiate risk level trajectory plans.

Discussion
This paper introduces a mathematical framework for driving risks assessment, and we demonstrate our work in typical conflict and risky scenarios on motorways. The approach can capture different lane markings, shoulder characteristics and variability of surrounding vehicle behaviours, and describe continuous risk for any 2D motion. Therefore, it can be generalised to highway driving.
We have relied on several simplifying assumptions to operationalise the framework. One should be cautious about the implications of these assumptions while using the approach and interpreting the results. The crash severity is expressed solely as the expected crash energy, an important characteristic of the crash process. But the outcome of a crash is influenced by several other factors such as vehicle's material stiffness at the point of impact, the presence of passive safety systems such as seat belts; and the material property of crash bumpers that are designed to partially dissipate the crash energy. Considering these elements in the crash severity quantification will improve the accuracy of the model. The crash between two vehicles is assumed to be inelastic, where the vehicles move together after the first contact. However, in reality, vehicle crashes involve some amount of elasticity, where vehicles rebound from each other. A scaling coefficient can be introduced to capture this effect. The acceleration variability distributions in the forward and lateral direction were assumed to be independent. However, they can be correlated during manoeuvres such as lane changes. Likewise the acceleration in different time steps was assumed to be uncorrelated. If such a correlation is detected while estimating the parameters, the variability should be modelled as a joint distribution.
One noteworthy issue about the model parameters is that the mean and especially the standard deviation of the normal distribution may influence on the prediction performance. In the cut-in scenario in simulation, a too-large standard deviation of acceleration may lead to false warnings. In the hard-brake scenario, a too-small standard deviation may result in false negative cases especially when the relative velocity of the two vehicles is negligible. In general, the model is more sensitive to the standard deviation than to the mean of acceleration. It is recommended to calibrate the mean and standard deviation with empirical data in different driving scenes and use motion cues such as lane changes to predict the evolution of mean and standard deviation evolution in time.
The crash probability estimation is based on the assumption that the subject vehicle maintains its current state; however, it might not do so during a system failure. The probability of system failure is a critical dimension of the risk measure. This dimension can be incorporated in the present framework treating the future motion of the subject vehicle as a stochastic feature that is related to the probability of a system failure. The probability of system failure is influenced by various factors such as hardware faults, software failures, which can be modelled as a fault tree (Wang et al., 2010). The accuracy of crash probability estimate can be improved by accounting for vehicle yawing in the motion prediction algorithm.
The additivity property of the artificial field is constructed under the condition that driving risk with two distinct road objects is independent. This condition is not met in a conflict involving multiple vehicles when the mutual vehicle interactions could influence each other's motion. However, such events are less prevalent than crashes between a vehicle pair.
For illustrative purpose, the risk measures in this paper were generated using parameter values selected from the literature. The accuracy of the risk measure can be improved by estimating them directly from a trajectory dataset or via test experiments as follows: • The parameter k, in potential PDRF, can be estimated in multiple ways. An objective approach is based on the relationship between k and the coefficient of restitution e res as follows: = k e (1 ) res 2 . The coefficient of restitution can be measured from vehicle crash tests with the specific road object (Noon, 1994). Some studies report that e res increases exponentially with the impact velocity (Noon, 1994). However, identifying e res for multiple combinations of the vehicle (at different velocity ranges) and boundary objects that appear in the dataset would be practically strenuous. Another, rather subjective, approach is to use the odds of injury (Zou et al., 2014) as the value of k. Such an approach can loosely capture the influence of the energy dissipation capacity of boundary objects on crash severity.
• The parameters of acceleration variability distributions can be estimated from the acceleration samples of a vehicle during a finite duration according to Eq. (5). The acceleration samples of a vehicle can be extracted from its trajectory data (simulated or video based naturalistic data) or can be acquired directly from the vehicle's data log. The sampling horizon should be small enough to capture the local traffic conditions and long enough to ensure statistical validity. Further research is needed to find an optimal sampling horizon. The acceleration noise when a vehicle is stopped in traffic is zero, which might distort its estimated value. Therefore, the acceleration noise should be measured only while the vehicle is moving. Furthermore, manoeuvres such as lane changing, turning and evasive braking are expected to feature distinctive sets of parameters, which require further study.  F.A. Mullakkal-Babu, et al. Transportation Research Part C 118 (2020) 102716

Conclusion
We presented an approach to assess driving risk, which employs a probabilistic motion prediction scheme within the framework of artificial potential field theory. The approach was designed to yield a continuous risk estimate and to account for important aspects of the driving risk: the crash severity and motion uncertainty. As a proof of concept, we illustrated the application of the approach with examples of safety assessment problems. Firstly, we applied the approach to analyse the risk of three near-crash situations selected from a naturalistic dataset. Here, we employed the risk estimate that describes the risk of a crash at a finite future time step. It was observed that the risk description qualitatively reflects the document narration of the situation; the evolution of risk in the situation. Secondly, we evaluated the accuracy of single-step PDRF ( = 3 s) in identifying the risk in two safety-critical situations, due to a cut-in neighbor and hard braking preceding vehicle. We observed that PDRF risk estimate was more accurate than TTC in both scenarios. Finally, we applied the model to estimate the risk of four possible trajectories that the subject could pursue while approaching a typical lane drop section. Here we employed the multi-step risk estimate that describes the risk of a crash at a sequence of multiple future time steps. It was observed that the risk estimates of the trajectory plans are plausible and clearly mark the safest trajectory, which is consistent with the well known Surrogate Measure of Safety: generalised TTC.
Both sets of illustrations demonstrated certain properties of the driving risk estimate. The risk descriptions vary within a finite intuitive range, i.e. between 0 and the expected crash energy. The results of both example sets showed that the risk trends described by the PDRF model, in general, were consistent with the prominent SMoS: Time To Collision. However, the risk measure (strength of PDRF) contained additional information: crash severity, which was seen to be the factor that differentiates the risk levels of trajectories. Based on the exhibited properties, and the relevance of the component aspects in characterising the driving risk, the proposed approach can be applied in analysing the safety of vehicle interactions and as a risk estimate in path planning algorithms (Mullakkal-Babu et al., 2020). As part of the assessment approach, we presented a probabilistic motion prediction scheme, which employs a distribution of acceleration variation to approximate the motion uncertainty. The use of acceleration variability makes the prediction scheme analytically tractable. The parameters of the distribution could be measured by monitoring the vehicle for the finite duration; and they are known to be sensitive to factors relevant to safety analysis such as road geometry, driver aggressiveness and traffic congestion.
The examples presented in this work are illustrative. The single step approach should be tested with large-scale naturalistic data (including crash, near-miss and regular driving) to evaluate its effectiveness in detecting risky situations. Such a test is also necessary to identify the optimal value of that provides the highest number of accurate detection. Similarly, the effectiveness of the multi-step approach as an ex-ante risk estimate should be evaluated by employing it in an advanced vehicle control design. We also acknowledge that the current model does not capture vulnerable users, as the behavioural uncertainties of these users are significantly different from that of vehicles and deserve separate research attention. We leave this as one of our future research directions.