Behavioral event detection and rate estimation for autonomous vehicle evaluation

An “autonomous driving system" (ADS) refers to a level 3, 4, or 5 system [12] that is capable of performing a dynamic driving task on a sustained basis. It is comprised of both software and hardware components, and any vehicle equipped with this system may be referred to as an "Autonomous Vehicle" (AV).

Development of safe and high quality ADSs is a very active area of research within both academia and industry [15, 22]. Waymo is a prominent player in this space, with its origins in 2009 as the Google Self-Driving Car Project. Today Waymo produces Level 4 ADSs which currently provide fully autonomous rides, i.e., with no human in the driver seat, to members of the public in Arizona and in California. Waymo also has parallel work in the autonomous trucking and delivery space.

At Waymo, rigorous validation and evaluation of the ADS performance is a core tenet of the deployment philosophy. Prior to releasing AVs equipped with a new ADS onto public roads, we must have confidence that the safety and other performance metrics meet our stringent quality bar. Previous Waymo publications have already outlined the safety and readiness determinations [20, 13]. Complementing those publications, here we articulate some of the statistical and data science challenges that arise during the evaluation process.

ADS performance evaluation typically involves three primary components:

1. (1)

   Definition of specific undesirable¹¹1There are many desirable behavioral events that are also of interest to track. For simplicity, in this paper we assume that these events can all be negated to produce a corresponding undesirable event. E.g., a desirable event may be “sufficiently pullover to pick up passenger out of the travel lane,” and the associated negative metrics would be things like “pullover had X percent overlap with travel lane” or “failed to pullover.” behavioral events, along with associated severity levels or other relevant annotations.

2. (2)

   Rate estimation for the behavioral events’ relative frequencies on-road, typically in terms of the number of incidents per driving mile (IPM).

3. (3)

   Accommodation of product requirements from the stakeholders who use these data for decision-making.

More details on each of these components are provided below.

The behavioral events of interest referenced in component (1) span a spectrum from critical safety situations such as collisions (very rare) to behavioral quirks such as unnecessary slowing (more common). While on its surface this aspect of performance evaluation appears straightforward, in practice rigorous definitions and delineation of edge cases can be a challenge. These challenges are left as out of scope for the current paper. Rather, we focus on data and statistical questions associated with the components (2) and (3), and how we are addressing these challenges at Waymo.

There are several data and statistical questions that arise when producing the rate estimates referenced in component (2):

1. a)

   What data is best suited to produce these rate estimates? Waymo has autonomously driven tens of millions of miles on public roads as well as tens of billions of miles in simulation. The former have minimal observation error but may have sparse coverage of rare events; the latter may be subject to simulation artifacts but can more easily ensure coverage of rare events.

2. b)

   Given the driving data, how can we reliably count these behavioral events? The ease of measurement for an individual class depends not only on the event’s frequency (denser events are easier to measure empirically on smaller datasets), but also our ability to algorithmically identify relevant events in the driving data. Identification of behavioral events in a purely algorithmic fashion is difficult and usually requires human review, often of varying levels of expertise ("tiered event review").

3. c)

   How can we quantify statistical uncertainty via reliable confidence intervals for the rate estimates? IPM is not a typical metric of interest within the technological industry outside of AVs, nor is it a quantity that receives broad attention from the academic literature. Moreover, the uncertainty introduced from imperfect event detection in (b) further complicates the statistical theory. Finally, these confidence intervals must accommodate a broad range of events, including rare event scenarios, for which traditional asymptotic methods can fail easily.

These statistical challenges are further constrained by the product requirements referred to in component (3):

• •

  In order to drive efficient decision-making at Waymo, these rate estimates must be available as soon as possible after mileage collection. As such, we need the ability to compute rate estimates even when not all candidate events have been reviewed. Early rate estimates may have greater uncertainty than the final estimates, but they should still be unbiased and provide robust confidence intervals.

• •

  The most practically important event classes are very sparse (< 5-10 events), which means the most challenging statistical scenarios are also the most important.

The rest of this paper outlines how we address 2(b) and 2(c) while accommodating these product constraints. In Section 2 we introduce a streaming partial tiered event review framework. Then, Sections 3, 4, 5 outline the data generation model, estimation, and confidence interval construction implied by this framework. Section 6 provides numerical studies to assess the performance of the proposed methods for both rare and common event cases. Finally, Section 7 concludes the paper with some discussion. Technical details are provided in the appendix, including a relation between Poisson, multivariate Hypergeometric and multinomial distributions, which may be of independent interest.

As described above, computing the rate at which a behavioral event occurs will require us to first identify and count those events within our driving data. To facilitate this detection at scale, we rely on a set of sophisticated algorithms that can automatically identify events of interest. These classification problems are conceptually similar to many other data science applications in industry. However, unlike most other application settings, edge cases are incredibly important within the autonomous vehicles industry and cannot be swept aside or averaged out. These edge cases, for us, are frequently the most interesting and most important events.

Take, for example, safety-relevant behavioral events which occur very rarely. Errors in their classification will show up as a very small percentage of the algorithm’s aggregate performance metrics, but lack of detection for these events can result in a very serious lack of information when deciding whether to put an ADS on public roads. This risk is compounded by the fact that estimating the recall loss by these algorithms can be very challenging and time-consuming since it requires finding sparse False Negatives in a very large collection of True Negatives.

To ensure confidence in our safety estimates at Waymo, we supplement these algorithmic solutions with a final stage of human review. When tuning any algorithmic classifier, one can adjust the operation point to favor either precision or recall, with very few algorithms being able to achieve high performance in both simultaneously. Our solution allows us to tune the algorithms for very high recall, and then lean on a final layer of human reviewers to fix any gaps in precision. Combined, the algorithm and human reviewers represent a multi-stage strategy that is both scalable and provides trustworthy rate estimates.

Let $\mathcal{E}_{0}$ refer to the corpus of candidate events identified by the algorithms. For simplicity, we assume the algorithm feeding events to the human reviewers has perfect recall, so the human reviewers only need to address precision. Then, for any candidate event in $\mathcal{E}_{0}$, we ultimately want to confirm whether the event was a True Positive or False Positive detection of the behavioral event of interest.²²2Note that detailed instructions and rigorous training are often needed to ensure that human reviewers provide high quality labels, particularly for the most nuanced behavioral events. For the purposes of this paper, we assume the human labels are perfect and, as such, recall from this system is perfect.

A well-optimized event review configuration should maintain both high precision and high recall for True Positive events. If the precision is low, then the rate estimates will be too high; if the recall is low, then the rate estimates will be too low. With this as our motivation, we have developed an event review configuration at Waymo referred to as streaming partial tiered event review. The core concepts motivating this configuration are outlined below.

Under a complete event review configuration, all candidate events would be reviewed. However, for many applications this can be inefficient and labor intensive. Instead, reviewing a subset of the candidate events, i.e., partial event review, is often sufficient to estimate event rates with acceptable levels of uncertainty. Beyond efficiency, supporting estimates based on partial event review also ensures that we can fulfill the product requirements by providing trustworthy estimates at any stage of the event review process (e.g., checking results after only 10% of events have received review), and that the uncertainty in those estimates will continually be reduced as more events are reviewed.³³3In addition to ongoing event review, the on-road ADS logs and simulated autonomous driving provide a continuous stream of candidate events that are being fed into the event review pipeline. As a result, the corpus of candidate events is constantly changing. For simplicity, this is not addressed here, and we instead assume the full corpus of candidate events is available and unchanging prior to the event review initiation.

In its simplest form, partial event review utilizes uniform random sampling to select a subset of candidate events for labeling, but this is unlikely to be the most efficient strategy for most use-cases. To reduce sampling error, stratified random sampling can instead be used; this involves partitioning candidate events into mutually exclusive, cumulatively exhaustive strata based on pertinent event metadata (e.g. road characteristics, algorithmic model scores). The sampling probability may be identical across strata (i.e. proportional allocation), or may differ across strata (e.g. Neyman allocation or alternative allocation schemes). In some sections of this paper we simplify to assume a single stratum, without loss of generality, while in other sections we outline the statistical implications of this stratification in the data.

For operational efficiency, an event review configuration may utilize several tiers of human reviewers with increasing levels of expertise in the higher tiers. This tiered configuration is arranged such that each tier improves the precision of the candidate event corpus while maintaining recall. For example, in a three tiered process, the first tier would exclude clear false positive candidate events while escalating all other events to the second tier. The second tier would exclude additional false positive candidate events, while escalating the remainder to the third tier. The third tier would then make the final determination of false positive vs. true positive status. Each tier may itself consist of evaluation by a single reviewer or multiple reviewers who reach consensus.

Operationally, partial tiered event review can be implemented as either a batched or streaming process. Under a batched process, any candidate event that is reviewed in Tier 1 will complete its journey through the system and be confirmed as either a False Positive or make it to the final tier and be confirmed as a True Positive. In this setting, the rate estimates do not need to account for the details of how candidate events are escalated between tiers; the rate estimates are ultimately identical between a partial tiered setting and a partial untiered setting. However, at Waymo it is often the case that we need to produce unbiased rate estimates while many events are still only partway through the event review pipeline. E.g., a candidate event may have already been reviewed by Tier 1 and escalated, but hasn’t received further review yet. This is considered to be a streaming process where rate estimates are made based on the data available at any snapshot in time.

Let $e_{0}=||\mathcal{E}_{0}||$ be the number of candidate events in the full corpus.

Now suppose we sample from the candidate events via a stratified sampling scheme with $H$ comprehensive and mutually exclusive strata. We can generalize all the notations in Section 3.1 in a straightforward way by an additional subscript $h\in\{1\cdots H\}$ indicating the stratum. The detailed data generation model is provided as Algorithm 2, which calls Algorithm 1 for each stratum $h$ by taking $\{\lambda_{ht}:0\leq t\leq T\}$ (associated with latent variables) and $\{\pi_{ht}:1\leq t\leq T\}$ (associated with tier level sampling) as the input parameters, and outputing the tier level sample sizes ($n_{ht}$) and the numbers of tier level escalations ($e_{ht}$) as the observed data. The parameter of interest can be written as

With the mileage $m$ fixed, conceptually, the model for the data generation can be written as $p(\mathcal{D}|\{m,H,T,\{\lambda_{ht},\pi_{ht}\}\})$, where $\mathcal{D}$ consists of $\{e_{ht}:1\leq h\leq H,0\leq t\leq T\}$ and $\{n_{ht}:1\leq h\leq H,1\leq t\leq T\}$. The detailed data generation procedure is provided by Algorithm 2. This is a parametric model, and the parametric bootstrap CI [6] can be constructed as follows:

Let $\{\hat{\lambda}_{ht},\hat{\pi}_{ht}\}$ be the point estimates of $\{\lambda_{ht},\pi_{ht}\}$ respectively, where $\hat{\lambda}_{ht}$s are obtained by applying the MLE estimator described in Theorem 1 to each stratum $h$ separately, and the following empirical estimate is used for $\pi_{ht}$ according to (3.2):

which is also MLE when $n_{ht}>1$. It may be worth noting that when $n_{ht}=1$ and $e_{h(t-1)}>1$, the MLE for $\pi_{ht}$ by the sampling model (3.2) is simply 0, which is close to (5.1) for large $e_{h(t-1)}$ but would be less desirable in general.

Then draw $B$ i.i.d. samples of $\mathcal{D}$ according to the generative process outlined in Algorithm 2 by plugging in the parameter estimates from the observed data. For $b=1,\cdots,B$,

For each draw, calculate the parameter estimates $\{\hat{\lambda}_{ht}^{(b)},\hat{\pi}_{ht}^{(b)}\}$ by replacing $\{e_{ht},n_{ht}\}$ in (4.3) and (5.1) with corresponding quantities in $\mathcal{D}^{(b)}$, and let

The $1-\alpha$ confidence interval for $\theta$ is given by the $\alpha/2$ and $1-\alpha/2$ quantiles of $\{\hat{\theta}^{(b)}:b=1,\cdots,B\}$.

In order to construct the Wald CI, we need to estimate $var(\hat{\theta})$. Unfortunately, the exact formula for $var(\hat{\theta})$ is extremely complicated. Instead, we derive a first order approximation by assuming $m$ is large.

It may be interesting to note that $\hat{\lambda}_{hT}$ as defined by (4.1) with index $h$ added to indicate the stratum can be rewritten as

where

and $\hat{\pi}_{ht}$s are defined by (5.1). Thus $\hat{\theta}$ in (4.4) can be rewritten as

For large $m$, we may approximate the model (3.2) with subscript $h$ added to indicate stratum as

Then by iteratively applying Lemma 2, we have

where $\pi_{h}=\prod_{t=1}^{T}\pi_{ht}$. By combining this with (5.2), we can obtain the following result.

As $m\rightarrow\infty,$ $\hat{\pi}_{ht}\rightarrow\pi_{ht}$ and thus $\hat{\pi}_{h}\rightarrow\pi_{h}$, the proof follows from Slutsky’s theorem.

From Theorem 3, the Wald CI can then be constructed by $\hat{\theta}\pm z_{\alpha/2}\sqrt{\frac{1}{m}\sum_{h=1}^{H}\hat{\lambda}_{hT}\hat{\pi}_{h}^{-1}}$ based on the plug-in principle, where $z_{\alpha/2}$ is the $1-\alpha/2$ quantile of $\mathcal{N}(0,1)$.

As is well known, proper coverage from the Wald CI requires the asymptotic normality to hold well. As such, it may not be desirable when the expected event counts are low (i.e., for rare events or for low mileage); for example, its lower bound may become negative.

Given the limitations of the normal asymptotic approximation outlined above, we propose an alternative approach to estimating the CIs for $\hat{\theta}$. Consider the estimate in (5.2), which can be approximated by the quantity below:

where by (5.3) the right hand side is a weighted sum of independent Poissons (WSIP). The approximate WSIP becomes clear by rewriting $\hat{\theta}$ as,

The distribution for WSIP random variables has been studied extensively in the statistics literature [10, 16, 17, 8]. A WSIP is close to normally distributed when there are high event counts, but can be far from normal for rare events. This again motivates the consideration of a CI methodology that does not rely on normality. In particular, the original Gamma method proposed by [7] has been conjectured and also shown numerically to provide good coverage properties for the distribution mean in all simulated scenarios [10, 8].

Our third CI methodology leverages these WSIP properties from the literature and can be described as follows:

• •

  The lower bound is the $\alpha/2$ quantile of the Gamma distribution with mean $\hat{\theta}$ and variance $\sum_{h=1}^{H}w_{h}^{2}e_{hT}$.

• •

  The upper bound is the $1-\alpha/2$ quantile of the Gamma distribution with mean $\hat{\theta}+w_{M}$ and variance $\sum_{h=1}^{H}w_{h}^{2}e_{hT}+w_{M}^{2}$, where $w_{M}=\max_{1\leq h\leq H}w_{h}$.

Note that the above is based on the original Gamma method, but one may also apply other Gamma methods for WSIP such as the mid-p Gamma interval [8].

We provide several simulation studies to assess and compare the performance of the three confidence interval methods defined in Section 5. The first two simulation studies compare the performance when the behavioral event is rare vs. when the behavioral event is common, using fixed values for most parameters. The third study illustrates that these coverage properties hold for a wide range of parameter values, beyond the specific values chosen for the first two examples. Throughout the simulation section we set $m=1$ for convenience (the unit could be 1 mile, 10K miles, or one million miles).

In each study, we simulate the data in a consistent manner. Following Algorithm 2, we generate the $x_{hs}\sim Poisson(\lambda_{hs}),h=1,\ldots,H$ and $s=0,\dots,T$. We set $H=5$ strata and assume $T=3$ tiers. Then $e_{h0}=\sum_{s=0}^{T}x_{hs}$ will be the total number of candidate events entering the tiered event review pipeline in stratum $h$.

For the first two simulation studies in Section 6.2 and 6.3, we set the sampling rates

Note that when this event review methodology is put into practice one can imagine a fixed set of events to review, and the Tier 1 sampling rate will steadily increase over time as reviewers complete more events. We simulate this event review progression by varying the sampling rate at Tier 1, $\pi_{1}\in[0.1,1.0]$. In this way, we can verify whether the confidence intervals have good coverage properties throughout the entire event review process.

For the third simulation study in Section 6.4, we consider a wide range of scenarios. For each scenario, we generate each Tier 1 sampling rate $\pi_{h1}\sim\text{Uniform }(0,1)$, and generate $\pi_{ht}\sim\text{Uniform }(\pi_{h(t-1)},1),t=2,3,$ so that the sampling rate for each strata will not decrease as Tier goes up. $\{\lambda_{ht}\}$ is specified later.

Given each set of parameter values ($\{\lambda_{ht}\}$ and $\{\pi_{ht}\}$), we simulate the tiered partial event review data $\{e_{ht}\}$ and $\{n_{ht}\}$ following Algorithm 2 and replicate 1000 times. This allows us to obtain the empirical coverage and average CI width under each set of parameters.

Some behavioral events are relatively common, and the observed event counts will be relatively high. As described earlier, this can result in better behaved quantities that more closely follow a normal distribution. Here we assess the performance of our proposed CI methodologies in this setting.

To represent this common event scenario, the values below are used for $\{\lambda_{ht}\}$:

With this parameter setting, $\theta=58$.

Figure 6.1 illustrates the coverage properties for the three CI methods when the behavioral event is common. The x-axis represents the sampling rate, allowing us to assess how well these CIs fulfill the product requirement of adequate coverage regardless of the sampling rate. For all three methods, the coverage rates are more biased (over- or under-covering) when sampling rates are low, with improvements in coverage as event review progresses. As sampling rates increase, the parametric bootstrap rapidly converges close to expected coverage rates, and also produces the narrowest intervals. However, the under-coverage at low sampling rates is non-negligible (82% coverage vs. the desired 90%). In contrast, the Gamma CI consistently over-covers (98% coverage for low sampling rates) and produces the widest intervals. However, the Gamma method also provides the most symmetric errors, with similar error rates for the upper and lower bounds; this guarantees that the intervals from this approach will not systematically under- or over-estimate the rates. Finally, the Wald CI falls between the other two methods, with intervals that somewhat over-cover.

Next, we consider more rare events, where the observed event counts will be quite low. In this setting, the estimates are not going to cleanly follow a normal distribution, making the CI calculation potentially more challenging. Here we assess the performance of our proposed CI methodologies in this setting, again tracking how the coverage properties evolve as more events receive review.

To represent this scenario, we use the same $\lambda_{ht}$ as before for $t=0,1,2$ but reduce $\lambda_{h3}$ for all $h$:

With this parameter setting, $\theta=11$.

Figure 6.2 reports the coverage properties for the three confidence interval methods when the behavioral event is rare. The patterns in this figure are largely similar to the patterns in Figure 6.1. The primary difference between the rare event case and the common event case is the performance at low sampling rates. The under-coverage by the Wald and Parametric Bootstrap is severe at the lowest sampling rates, off by as much as 30 percentage points (60% coverage vs. the desired 90%). In this case, the over-coverage produced by the Gamma CI (97% vs. the desired 90%) is clearly preferable. That said, the Gamma CI is much wider than the Parametric Bootstrap, as much as 3x wider for low sampling rates.

In each of the simulations above, the majority of parameters were fixed, while the first tier sampling rate was allowed to vary. In our final simulation we continue to vary the sampling rate, as described earlier, but consider $10^{5}$ different scenarios covering a wide range of parameter values $\lambda_{ht}$. In each scenario, we simulate the underlying Poisson rates as follows:

where $\mu_{ht}\sim\text{Uniform}(1,4)$ independently. This simulation study is inspired by the simulation setup in [8].

Figure 6.3 reports the performance comparison. Unlike Figures 6.1 and 6.2 where the x-axis represented a range of sampling rates, the x-axis here corresponds to the expected observed number of true positive events per mile, i.e. $\frac{1}{m}\sum_{h}E(e_{hT})$ which equals $\sum_{h}\lambda_{hT}\prod_{t=1}^{T}\pi_{ht}$. In the earlier examples, the underlying event rates were fixed and only the sampling rates varied; in this study, both the sampling rates and the event rates are varying simultaneously. In this setup, this x-axis allows for consideration of how the coverage properties evolve in the asymptotic sense, regardless of whether it is increased review or increased event occurrence.

The pattern is overall similar to the previous examples, for example, for the lower bound, both Gamma method and Wald method maintain the right coverage. However, the Wald method can be badly under-covered for its upper bound, while parametric bootstrap can be badly under-covered for both lower bound and upper bound. On the other hand, the Gamma method is the only method that always produces reliable confidence intervals for all scenarios. Not surprisingly, this method also produces the widest intervals on average.

Consider an urn of balls, where each ball has a color $k$ with $k\in\{0,1,\cdots,K\}$. Let $z_{k}$ be the number of balls with color $k$, and let $s_{z}=\sum_{k=0}^{K}z_{k}$ and $\vec{z}$ be the column vector of $(z_{0},\cdots,z_{K})$.

Let $n\leq s_{z}$ be a random number whose distribution only depends on $s_{z}$.

Randomly draw $n\geq 1$ balls from the urn without replacement, and let $y_{k}$ be the number of balls with color $k$ for $k\in\{1,\cdots,K\}$. Let $s_{y}=\sum_{k=1}^{K}y_{k}$ and $\vec{y}$ be the column vector of $(y_{1},\cdots,y_{K})$.

Then, conditional on $(\vec{z},n)$, $\vec{y}$ follows a multivariate Hypergeometric distribution:

The data for a single stratum can be described as $(e_{0},n_{1},e_{1},\cdots,n_{T},e_{T})$. Note that in case the partial event review process is terminated earlier due to no escalation, say $e_{t}=0$ for some tier $t<T$, then there will be no data from higher tiers (i.e., $n_{s}=0$ and $e_{s}=0$ for $s>t$ according to Algorithm 1), so we may simply reset $T$ to $t$ for the data likelihood calculation. Let $\vec{e}=(e_{0},\cdots,e_{T})$ and $\vec{n}=(n_{1},\cdots,n_{T})$. The observed data consists of $\vec{e}$ and $\vec{n}$. Without loss of generality, we will take the mileage $m=1$.

To follow Algorithm 1, let $\mathcal{E}_{0}$ consist of all $e_{0}$ candidate events. Let $\mathcal{E}_{t}$ be the set of $e_{t}$ events which are escalated by Tier-$t$ for $t=1,\cdots,T$, so $\mathcal{E}_{T}$ consists of the $e_{T}$ events which are True Positive. Furthermore, each candidate event is associated with a label, which is either FP-$t$ (if the event would be escalated by Tier-$(t-1)$ but rejected as False Positive by Tier-$t$) for some $t\in\{1,\cdots,T\}$, or TP (if the event would not be rejected by any Tier), under the complete event review configuration.

For $s\in\{0,\cdots,T-1\}$, let $x_{ts}$ be the number of FP-$t$ events in $\mathcal{E}_{s}$, $\forall s\leq t<T$, and let $x_{Ts}$ be the number of TP events in $\mathcal{E}_{s}$, then

For convenience, we also introduce $x_{TT}\equiv e_{T}$. With some abuse of notation, $(x_{00},\cdots,x_{T0})$ is identical to $(x_{0},\cdots,x_{T})$ as defined by (3.1).

Let $\vec{x}_{s}=(x_{ss},\cdots,x_{Ts})$, for $s=0,\cdots,T$. Table 1 summarizes the notation of $\{x_{ts}:0\leq s\leq t\leq T\}$ and relation with $\vec{e}$.

According to the sampling procedure described in Algorithm 1, the following conditional independence holds:

Also note that $n_{s}$ only depends on $e_{s-1}$. Therefore the full likelihood function of $\{\vec{x}_{s}:0\leq s\leq T\}$ and $\{n_{s}:1\leq s\leq T\}$ can be written as

According to the Poisson assumption described in Algorithm 1,

According to the sampling procedure,