Development of Analysis Methods Using Recent Data (2012)

Below is the uncorrected machine-read text of this chapter, intended to provide our own search engines and external engines with highly rich, chapter-representative searchable text of each book. Because it is UNCORRECTED material, please consider the following text as a useful but insufficient proxy for the authoritative book pages.

10 C H A P T E R 2 State-Space Model Chapter 1 illustrated the trajectory-based approach with simple braking-to-stop models estimated from data extracted from site-based video. An important objective of this proj- ect was to develop a common analytic framework that could be applied to data from either site- or vehicle-based sensor configurations and could accommodate events more com- plicated than constant speed followed by simple braking- to-stop. Ultimately, the goal is to model vehicle trajectories in two dimensions, where both braking and steering could serve as evasive actions. This report, however, focuses on the simpler problem of modeling vehicle motion in one direction. The basic idea is to model driver behavior as a piecewise constant series of acceleration changes, which are then treated as inputs into a dynamic trajectory model. The vehicle’s state at a given time is its location and velocity, and the trajectory model takes the acceleration input sequence and numeri- cally integrates the associated differential or difference equations to produce time histories of vehicle locations and speeds. For discrete-time data, the trajectory model can be conveniently represented using the generic linear state-space form x t Ax t Ba t y t Cx t +( ) = ( )+ ( ) ( ) = ( ) 1 2 1( . ) where x(t) is a vector of state variables (position and velocity), a(t) is a vector of input variables (accelerations), and y(t) is the vector of observed variables. A, B, and C stand for matrices of coefficients. The nature of A, B, C, x(t), and y(t) will vary depending on the class of events being modeled and the sort of data avail- able. For two vehicles following on a straight road, the sim- plest trajectory model consists of two state variables for each vehicle, location and speed, with linear acceleration values as inputs. That is, if D denotes the basic time interval of the data, then the deterministic progression for a leading and follow- ing vehicle can be captured by the linear equation x t x t v t v t 1 1 2 1 1 1 2 1 1 0+( ) +( ) +( ) +( )       = ∆ ∆       ( ) ( ) ( ) 0 0 1 0 0 0 1 0 0 0 0 1 1 2 1 • x t x t v t v t a t 2 0 0 0 0 0 0 1 ( )       + ∆ ∆       ( ) ( )    a t2 2 2( . ) Here, x1(t) and x2(t) give the locations of the leading and following vehicles at time t, v1(t) and v2(t) are the corre- sponding speeds, and a1(t) and a2(t) are the accelerations. For motion in two directions, a similar structure can be used but with state and input variables for each direction. Given initial values for the state variables and the time history for the inputs, the trajectories of both vehicles can be replicated. Since Equation 2.2 describes vehicle motion irrespective of the data-collection scenario, the primary difference between the site- and vehicle-based scenarios will be the observation equation. For vehicle trajectories extracted from video, the observations consist of measurements of position for each vehicle, leading to an observation equation of the form y t y t x t x t v 1 2 1 0 0 0 0 1 0 0 1 2( ) ( )    =     ( ) ( ) • 1 2 2 3 t v t ( ) ( )       ( . ) A standard calculation shows that this system is observ- able, and the two-vehicle system decomposes into two one-vehicle systems. For the vehicle-based data from the 100-car study, the pri- mary observations are the speed of the following vehicle obtained from its speedometer, and the range and range rate Overview of Analytic Method

11 for the lead vehicle relative to the follower, obtained from the forward radar. This leads to an observation equation of the form y t y t y t 1 2 3 0 0 1 0 1 1 0 0 0 0 1 1 ( ) ( ) ( )       = − −       ( ) ( ) ( ) ( )       • x t x t v t v t 1 2 1 2 2( . )4 This system is not observable. For instance, it is not pos- sible to obtain estimates of absolute position for each vehicle. It also does not decompose into two independent sub- systems, so that estimation and inference will generally require working with both vehicles. For site-based data from the CICAS radar system, observa- tions consist of both position and speed for individual vehi- cles, leading to an observation equation of the form y t y t y t y t 1 2 3 4 1 0 0 0 0 1 0 0 0 0 ( ) ( ) ( ) ( )       = 1 0 0 0 0 1 1 2 1 2       ( ) ( ) ( ) ( )    • x t x t v t v t     ( . )2 5 This also decomposes into two separate subsystems, one for each vehicle. Given estimates of a driver’s initial speed, the times at which he or she changed acceleration, and the correspond- ing accelerations, the differential equations can be solved to give predicted time histories of that vehicle’s position and speed and predicted values for the observations. Fitting a trajectory model then involves searching plausible combi- nations of values for these input quantities to find those that best account for the data. The counterfactual simulation needed to assess the degree to which a near crash might have been a crash, in turn, involves using probability distributions characterizing the residual uncertainty in the model param- eters as input to a Monte Carlo simulation, where the state equation is integrated using random draws from this distri- bution and the occurrence or nonoccurrence of a collision recorded. In the early stages of this study, several different approaches to implementing these steps were experimented with, including modeling the differential equation model using response-surface approximations (1), nonlinear least- squares estimation using asymptotic normal approxima- tions to characterize posterior uncertainty (2), and Bayesian analysis using Markov Chain Monte Carlo (MCMC) simula- tion. The WinBUGS software (3) can be used to implement MCMC estimation for a variety of relatively complicated models, but the project team’s initial attempts to implement this by directly coding the differential equations in WinBUGS led to estimation runs with excessively long time demands. This problem was circumvented by using the WinBUGS dif- ferential equation interface. This provides compiled proce- dures that can be included in a WinBUGS model specification, which numerically solves ordinary differential equations using Runge-Kutta methods. As is the standard practice in working the MCMC estimation, exploratory analyses were first conducted using frequentist methods, in this case non- linear least-squares, implemented using either MATLAB (4) or R (5). This was done to understand the complexity of the acceleration model suggested by a given data set and to obtain reasonable starting values for the MCMC simulation. Bayes estimates for model parameters were then computed using WinBUGS, and counterfactual simulation was carried out using the MCMC sample of the posterior distribution for these parameters. Illustrative Example The general approach is illustrated here using case 104119 from the 100-car study. This event was a potential near crash that involved a lead vehicle and a following vehicle successively braking to a stop, with the follower stopping short of collision. This was a vehicle-based study, with the following vehicle being the instrumented vehicle. The data employed in this analysis were the speedometer-measured speeds for the instrumented vehicle and range and range rate for the leading vehicle, obtained from the follower’s forward radar. Step 1: Graphical Inspection of Data Figure 2.1 shows the time history of the speedometer-measured speeds of the following vehicle converted to units of ft/s. The piecewise constant nature to this relation was characteristic of all speedometer data obtained for 100-car study cases. Inspection of Figure 2.1 suggested a two-phase model where the following driver was initially traveling at about 21 ft/s and accelerating until about 11 s from the start of the data series. He or she then decelerated at a roughly constant rate until coming to a stop. Figure 2.2 shows the time history of the range and range- rate data provided by the instrumented vehicle’s forward radar. The discontinuities in the time-series result from peri- ods of missing data, when the forward radar apparently lost the leading vehicle as a target. To get an initial sense of the leading driver’s actions, an approximate speed profile for the leading vehicle was often used, obtained by adding the range rate to the following vehicle’s speedometer data. A plot of these approximate speeds is shown in Figure 2.3. This

12 suggests a three-phase model where a period of constant acceleration is followed by a period of roughly constant speed, which in turn is followed by a period of constant decel- eration leading to a stop. Step 2: Nonlinear Least-Squares Estimation of Proposed Models It is necessary to estimate the following vehicle’s initial speed, its acceleration during the first phase, the time at which deceleration began, and the deceleration characteristic of the second phase. To estimate these, a MATLAB script was writ- ten that took trial values for these parameters as inputs and simulated the following vehicle’s position and speed over time by solving the differential equations using a simple Euler’s method. The difference between the simulated speeds and the speedometer speeds was computed for each time inter- val, and the squares of these differences were summed to produce a measure of fit between the model and the speed- ometer data. This script was then embedded in a numerical Figure 2.1. Speedometer speeds of following vehicle. -10 0 10 20 30 40 50 0 2 4 6 8 10 12 14 16 Time in sec feet -15 -10 -5 0 5 10 feet/sec range in feet range rate in feet/sec Note: Discontinuities are due to missing data. Figure 2.2. Range and range-rate data for the leading vehicle.

13 search procedure to find the parameter values that minimized the sum-of-squares. The resulting estimates and initial approximate standard errors are shown in Table 2.1. For the lead vehicle, the three-phase model was fit to the approximate speed data shown in Figure 2.3. The resulting nonlinear least-squares estimates and approximate standard errors are given in Table 2.2. Step 3: Bayes Estimation of Vehicle Models Final Bayes estimates were computed using the MCMC soft- ware WinBUGS. In essence, WinBUGS generates a simulated realization of a Markov chain whose stationary distribution is the same as the Bayesian posterior distribution of the model parameters given the data. For this, the data were the speed profile for the follower and the range and range-rate profiles for the leader. Model parameters consisted of accel- eration profiles and initial speeds for each vehicle, and the parameters for both vehicles were estimated in combination simultaneously. In the WinBUGS model, predicted values for the follower’s speed and the range and range rate for the Figure 2.3. Approximate speed of lead vehicle determined from follower speed and range rate. Table 2.1. Nonlinear Least-Squares Estimates for the Following Vehicle’s Acceleration Model (Case 104119) Parameter Estimate Standard Error Initial acceleration (ft/s2) 1.45 .02 Final acceleration (ft/s2) −9.9 .17 Transition time (seconds from start) 11.2 .04 Table 2.2. Nonlinear Least-Squares Estimates for Lead Vehicle Acceleration Model (Case 104119) Parameter Estimate Standard Error First acceleration (ft/s2) 2.36 .05 Second acceleration (ft/s2) −0.84 .23 Third acceleration (ft/s2) −9.86 .40 First change point (s) 6.75 .06 Second change point (s) 10.4 .09

14 Table 2.3. Posterior Summary for Trajectory Model Parameters Variable Mean Standard Deviation 2.5%ile 97.5%ile Following Vehicle Initial speed (ft/s) 20.6 0.24 20.1 21.0 First acceleration (ft/s2) 1.45 0.04 1.38 1.54 Second acceleration (ft/s2) −9.47 0.21 −9.89 −9.05 First change (s) 11.15 0.06 11.04 11.26 Leading Vehicle Initial speed (ft/s) 18.92 0.31 18.3 19.52 First acceleration (ft/s2) 3.34 0.14 3.08 3.62 Second acceleration (ft/s2) 0.62 0.05 0.51 0.72 Third acceleration (ft/s2) −11.94 0.32 −12.58 −11.34 First change (s) 3.47 0.13 3.22 3.73 Second change (s) 10.61 0.05 10.51 10.72 leader were computed by numerically solving the differential equations using ordinary differential equation interface. For this case, a 10,000-iteration burn-in period followed by a 70,000-iteration MCMC sample produced acceptable con- vergence. Table 2.3 summarizes the Bayes estimates for both leader and follower. Figure 2.4 shows the speedometer speeds of the following vehicle, together with the posterior means estimated from the MCMC sample. Figure 2.5 shows similar results for the range data from the forward radar. Both cases provide rea- sonable approximations to the observed information. Finally, Figure 2.6 shows the probability of a rear-end col- lision between the two vehicles as a function of counterfactual final decelerations by the following driver. These probabilities were computed using Balke and Pearl’s (6) Twin Network method, where the follower’s deceleration is set at a target value and then, for each outcome of the MCMC, sample val- ues for the remaining parameters are used as inputs to solve the differential equations. This describes what would have happened had the event involved those parameter values and the counterfactual follower’s deceleration. Simulated range values less than zero are taken to indicate a collision, and the fraction of the MCMC sample values that lead to collision is an estimate of the collision probability. In this case, had the follower braked at less than about –10 ft/s2, a collision would probably have resulted. 0 5 10 15 20 25 30 35 40 0 5 10 15 20 Time in sec Speed feet/sec Measured Predicted Figure 2.4. Measured following vehicle speeds and posterior mean predicted speeds from MCMC sample.

15 Extensibility. Statistics and Computing, Vol. 10, No. 4, 2000, pp. 325–337. 4. Using MATLAB (Version 6). The MathWorks, Inc., Natick, Mass., 2002. 5. R: A Language and Environment for Statistical Computing. R Foun- dation for Statistical Computing, Vienna, Austria, 2008. www .R-project.org. 6. Balke, A., and J. Pearl. Probabilistic Evaluation of Counterfactual Queries. Presented at 12th National Conference on Artificial Intel- ligence (AAAI-94), Seattle, Wash., 1994. References 1. Bayarri, M., J. Berger, D. Higdon, M. Kennedy, A. Kottas, R. Paulo, J. Sacks, J. Cafeo, J. Cavendish, C. H. Lin, and J. Tu. A Framework for Validation of Computer Models. Technical Report 128. National Institute of Statistical Sciences, Research Triangle Park, N.C., 2002. 2. Seber, G. A. F., and C. J. Wild. Nonlinear Regression, Wiley, New York, 1989. 3. Lunn, D., A. Thomas, N. Best, and D. Spiegelhalter. WinBUGS— A Bayesian Modelling Framework: Concepts, Structure, and 0 5 10 15 20 25 30 35 40 45 50 0 2 4 6 8 10 12 14 Time in sec Range in feet Measured Predicted Figure 2.5. Measured range and posterior mean predicted range from MCMC sample. 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 -32 -27 -22 -17 -12 -7 -2 Probability Deceleration (feet/sec2) Figure 2.6. Collision probability as a function of counterfactual values for follower’s final deceleration.