Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
148,768 characters · 34 sections · 104 citation commands
An Automatic Finite-Sample Robustness Metric: When Can Dropping a Little Data Make a Big Difference?
\notbool{arxiv}{
}
\seclabel{introduction} Ideally, policymakers will use economics research to inform decisions that affect people's livelihoods, health, and well-being. Yet study samples may differ from the target populations of these decisions in non-random ways, perhaps because of practical challenges in obtaining truly random samples, or because populations generally differ across time and place. When these deviations from the ideal random sampling exercise are small, one might think that the empirical conclusions would still hold in the populations affected by policy. It therefore seems prudent to ask whether a small percentage of a study's sample---or a handful of data points---has been instrumental in determining its findings. In this paper we provide a finite-sample, automatically-computable metric of how dropping a small amount of data can change empirical conclusions. We show that certain empirical results from high-profile studies in economics can be reversed by removing less than 1% of the sample even when standard errors are small, and we investigate why.
There are several reasons to care about whether empirical conclusions are substantially influenced by small percentages of the finite sample. In practice, even if we can sample from the population of direct interest, small percentages of the data are missing; either surveyors and implementers cannot find these individuals, or they refuse to answer our questions, or their answers get lost or garbled during data processing. As this missingness cannot safely be assumed random, researchers might care whether their substantive conclusions could conceivably be overturned by a missing handful of data points. Similarly, consumers of research who are concerned about potentially non-random errors in sample construction at any stage of the analysis might be interested in this metric as a measure of the exposure of a study's conclusions to this concern. Conclusions that are highly influenced by a small handful of data points are more exposed to adverse events or errors during data analysis, including p-hacking, even if these errors are unintentional.
Even if researchers could construct a perfectly random sample from a given study population, the target population for our policy decisions is almost always different from the study population, if only because the world may change in the time between the research and the decision. For this reason, social scientists often aspire to uncover generalizable or “externally valid” truths about the world and to make policy recommendations that would apply more broadly than to a single study population.
In this paper, we propose to directly measure the extent to which a small fraction of a data sample has influenced the central claims or conclusions of a study. For a particular fraction $\alpha$ (e.g., $\alpha = 0.001$), we propose to find the set of no more than $100 \alpha \%$ of all the observations that effects the greatest change in an estimator when those observations are removed from the sample, and to report this change. For example, suppose we were to find a statistically-significant average increase in household consumption after implementing some economic policy intervention. Further suppose that, by dropping 0.1% of the sample (often fewer than 10 data points), we instead find a statistically-significant average decrease in consumption. Then it would be challenging to argue that there is strong evidence that this intervention would yield consumption increases in even slightly different populations.
To quantify this sensitivity, one could consider every possible $1-\alpha$ fraction of the data, and re-run the original analysis on all of these data subsets. But this direct implementation is computationally prohibitive.\footnote{Indeed, young2019consistency finds it computationally prohibitive to re-run their analysis when leaving out every possible subset of two data points. To illustrate, consider an analysis that takes 1 second to run; checking removal of every 4 data points from a data set of size 400 would take over 33 years. See \secref{MIP} for more detail.} We propose a fast approximation that works for common estimators---including Generalized Methods of Moments (GMM), Ordinary Least Squares (OLS), Instrumental Variables (IV), Maximum Likelihood Estimators (MLE), Variational Bayes (VB), and all minimizers of smooth empirical loss (\secref{AMIP}). Computation of the approximation is fast, automatable, and easy to use, and we provide an R package on GitHub called “zaminfluence.”\footnote{\url{https://github.com/rgiordan/zaminfluence}. The name stands for “Z-estimator approximate maximum influence.” }
Our approximation is based on the classical “influence function,” which has been used many times in the literature to assess sensitivity to dropping one or a small number of datapoints (a discussion of related work can be found in \secref{related_work} below). However, prior work focused on outlier detection and visual diagnostics and considered small numbers of removed datapoints. In contrast, we relate the effect of ablating a non-vanishing proportion of datapoints to classical inference, with an interest in generalizing to unseen populations rather than detection of gross outliers, and analyze the accuracy of the empirical influence function as an approximation to leaving out a fixed proportion of data.
Specifically, we show that our approximation performs well using a combination of theoretical analyses, simulation studies, and applied examples. We demonstrate theoretically that the approximation error is low when the percentage of the sample removed is small (\secref{accuracy}). Moreover, for the cost of a single additional data analysis, we can provide an exact lower bound on the worst-case change in an analysis upon removing $100\alpha \%$ of the data (\secref{exact_lower_bound}). We check that our metric detects combinations of data points that reverse empirical conclusions when removed from real-life datasets (\secref{examples}). For example, in the Oregon Medicaid study finkelstein2012oregon, we can identify a subset containing less than 1% of the original data that controls the sign of the effects of Medicaid on certain health outcomes. In the Mexico microcredit study angelucci2015microcredit, we find a single observation, out of 16,500, that controls the sign of the ATE on household profit.
We investigate the source of this sensitivity when it arises, and we show that it is not captured in conventional standard errors. We find that a result's exposure to the influence of a small fraction of the sample need not reflect a model misspecification problem nor the presence of gross outliers. Sensitivity according to our metric can arise, even if the model is exactly correct and the data set arbitrarily large, if there is a low signal-to-noise ratio: that is, if the strength of the claim (signal) is small relative to a quantity that consistently estimates the standard deviation of the limiting distribution of root-$N$ times the quantity of interest (\secref{why}). For example, in OLS this “noise” is large when we have a high ratio of residual variance to regressor variance (\secref{influence_function_ols}). This noise can be large even when standard errors are small, because it does not disappear as $N$ grows. This result highlights the distinction between performing classical inference within a hypothetical perfect random resampling experiment, and attempting to generalize beyond the data to the world in which very small changes to the population are occurring over space and time.
We examine several applications from empirical economics papers and find that the sensitivity captured by our metric varies considerably across analyses in practice. In many cases, the sign and significance of certain estimated treatment effects can be reversed by dropping less than 1% of the sample, even when the t-statistics are very large and inference is very precise; see, e.g., the Oregon Medicaid RCT finkelstein2012oregon in \secref{example_medicaid}. In \secref{example_transfers}, we examine the Progresa Cash Transfers RCT angelucci2009indirect and show that trimming outliers in the outcome data does not necessarily reduce sensitivity. In \secref{example_microcredit_linear} we examine a simple two-parameter linear regression on seven Microcredit RCTs meager2020aggregating and, in \secref{example_microcredit_hierarchical}, we examine a Bayesian hierarchical analysis of the same data; these final two analyses show that neither very simple nor relatively complex Bayesian models are immune to sensitivity to dropping small fractions of the data. However, not all analyses we examine are non-robust. Certain results across the applications we examine are robust up to 5% and even 10% removal.
We recommend that researchers use our metric to complement standard errors and other robustness checks. Our goal is not to supplant other sensitivity analyses, but to provide an additional tool to be incorporated into a broader ecosystem of systematic stability analysis in data science yu:2013:stability. For example, since our approximation is fundamentally local due to the Taylor expansion, practitioners may also consider global sensitivity checks such as those proposed by leamer1984global, leamer1985sensitivity, sobol2001global,saltelli2004global, or the conventional breakdown frontiers approach of he1990tail, masten2020inference. Our method is also not a substitute for tailored robustness checks designed by researchers to investigate specific concerns about sensitivity of results to certain structures or assumptions. Applied researchers will always know more than econometricians about which specific threats to their empirical strategies are most worth investigating in order to solidify our trust in the results of any given analysis. And practitioners may well benefit from robustifying their analysis mostellertukeydata, hansen2008robustness,chen2011sensitivity even if they pass our check. Our metric is also complementary to classical gross error robustness (which we take to include outlier detection and breakdown point analyses) belsley:1980:regression,hampel1986robustbook. In particular, gross error sensitivity is designed to detect and accommodate arbitrary adversarial perturbations to the population distribution. We discuss similarities and differences between our work and other robustness measures in detail in \secref{related_work}.
We do not recommend researchers discard results that are not robust to removal of a small, highly-influential subset of data. While in certain cases such sensitivity may be concerning for specific, contextually-determined reasons, there is as yet no basis for doing so in general, as we have shown that such sensitivity can arise even if the conventional inference is valid in the strictest sense. However, we do suggest that researchers adjust their interpretation of results which are sensitive to dropping a small fraction of the data as being less generally applicable to somewhat differing populations, and less robust to minor corruptions of their random sampling assumption. Much as one would interpret statistically insignificant results as a failure to detect an effect rather than positively detecting the absence of an effect, sensitive results may indicate a failure to detect a transportable effect, but not necessarily a failure of classical inference in itself. We do not yet recommend any specific alterations to common inferential procedures based on our metric, but we believe this direction is promising for future research.
\seclabel{metric} \seclabel{MIP}
Suppose we observe $N$ data points $d_{1}, \ldots, d_{N}$. For instance, in a regression problem, the $n$-th data point might consist of covariates $x_n$ and response(s) $y_n$, with $d_n = (x_n,y_n)$. Consider a parameter $\theta \in \mathbb{R}^{P}$ of interest. Typically we estimate $\theta$ via some function $\hat{\theta}$ of our data. The central claim of an empirical economics paper is typically focused on some attribute of $\theta$, such as the sign or significance of a particular effect or quantity. A frequentist analyst might be worried if removing some small fraction $\alpha$ of the data were to
To capture these concerns, we define the following quantities:
In general, to compute the Maximum Influence Perturbation for some $\alpha$, we would need to enumerate every data subset that drops no more than 100$\alpha$% of the original data. And, for each such subset, we would need to re-run our entire data analysis. If $m$ is the greatest integer smaller than 100$\alpha$, then the number of such subsets is larger than $\binom{N}{m}$. For $N = 400$ and $m=4$, $\binom{N}{m} = 1.05 * 10^9$. So computing the Maximum Influence Perturbation in even this simple case requires re-running our data analysis over 1 billion times. If each data analysis took 1 second, computing the Maximum Influence Perturbation would take over 33 years to compute. Indeed, the Maximum Influence Perturbation, Most Influential Set, and Perturbation-Inducing Proportion may all be computationally prohibitive even for relatively small analyses.
To address this computational issue, we propose to instead use a (fast) approximation to the Maximum Influence Perturbation, Most Influential Set, and Perturbation-Inducing Proportion. We will see, for the cost of one additional data analysis, our approximation can provide a lower bound on the exact Maximum Influence Perturbation. More generally we provide theory and experiments to support the quality of our approximation. We provide open-source code\footnote{\url{https://github.com/rgiordano/zaminfluence}} and show that our approximation is fully automatable in practice (\secref{AMIP}).
We articulate our approximation in \secref{taylor_series,AMIP} below. First, in \secref{taylor_series} to follow we derive a Taylor series approximation to the act of leaving out datapoints. Though this approximation is based on a well-known first-order Taylor series approximation to the act of leaving out datapoints, known as the {\em empirical influence function} hampel1974influence,hampel1986robustbook, we will assume no familiarity with this work, deferring discussion of related literature to \secref{influence_function,related_work}. We then define our approximation to data dropping in \secref{AMIP}, using the observation that the finding the Maximum Influence Perturbation and its related quantities is trivial for the Taylor series approximation. We then conclude this section with some simple, concrete examples of our approximation in \secref{function_examples,linear_regression}.
\seclabel{taylor_series} We begin by a deriving a Taylor series approximation to the act of dropping data. Though this approximation is well-known as the empirical influence function (see \secref{influence_function} below for more details), we will derive the approximation assuming no prior knowledge other than ordinary multivariate calculus.
To form a Taylor series, we will naturally require certain aspects of our estimator to be differentiable. We now summarize common assumptions under which the Taylor expansion exists, and note that many common analyses satisfy these assumptions---including, but not limited to, typical settings for OLS, IV, GMM, MLE, and variational Bayes. Below, in \secref{accuracy}, we will state stricter sufficient conditions that guarantee not only the existence but also the finite-sample accuracy of our approximation.
For instance, the function that picks out the $p$-th effect from the vector $\theta$, $\phi(\theta) = \theta_{p}$, satisfies this assumption.
To form a Taylor series approximation to the act of leaving out datapoints, we introduce a vector of data weights, $\vec{w} = (w_1, \ldots, w_N)$, where $w_n$ is the weight for the $n$-th data point. We recover the original data set by giving every data point a weight of 1: $\vec{w} = \vec{1} = (1, \ldots, 1)$. We can denote a subset of the original data as follows: start with $\vec{w} = \vec{1}$; then, if the data point indexed by $n$ is left out, set $w_n = 0$. We can collect weightings corresponding to all data subsets that drop no more than 100$\alpha$% of the original data as follows:
Our approximation will be to form a Taylor expansion of our quantity of interest $\phi$ as a function of the weights, rather than recalculate $\phi$ for each data subset (i.e., for each reweighting).
To that end, we first reformulate our setup, now with the weights $\vec{w}$; note that we recover the original problem (for the full data) above by setting $\vec{w}=\vec{1}$ in what follows. Let $\hat{\theta}(\vec{w})$ be our parameter estimate at the weighted data set described by $\vec{w}$. Namely, $\hat{\theta}(\vec{w})$ is the solution to the weighted estimating equation
We allow that the quantity of interest $\phi$ may depend on $\vec{w}$ not only via the estimator $\theta$, so we optionally write $\phi(\theta, \vec{w})$ with $\phi(\cdot,\cdot): \mathbb{R}^{P} \times \mathbb{R}^N \rightarrow \mathbb{R}$. Whenever we write $\phi(\cdot)$ as a function of a single argument, we will implicitly mean $\phi(\cdot, \vec{1})$. We require that $\phi(\cdot,\cdot)$ be continuously differentiable in both its arguments. For instance, we can use $\phi(\theta,\vec{w}) = \theta_{p}$ to pick out the $p$-th component of $\theta$. Or, to consider questions of statistical significance, we may choose $\phi(\theta,\vec{w}) = \theta_{p} + 1.96 \sigma_{p}(\theta,\vec{w})$, where $\sigma_{p}(\theta,\vec{w})$ is an estimate of the standard error depending smoothly on $\theta$ and $\vec{w}$; this example is our motivation for allowing the more general $\vec{w}$ dependence in $\phi(\theta, \vec{w})$.
With this notation in hand, we can restate our original goal of computing the Most Influential Set as solving
Here we focus on positive changes in $\phi$ since negative changes can be found by reversing the sign of $\phi$ and using $-\phi$ instead. In particular, the zero indices of $\vec{w}^{**}$ correspond to the Most Influential Set: $\mis{\alpha} := \left\{n: \vec{w}^{**}_n = 0 \right\}$. And $\mip{\alpha} = \phi(\vec{w}^{**}) - \hat{\phi}$ is the Maximum Influence Perturbation. The Perturbation Inducing Proportion is the smallest $\alpha$ that induces a change of at least size $\Delta$: $\loprop{\Delta} := \inf\{ \alpha: \mip{\alpha} > \Delta\}$.
\seclabel{AMIP} Our approximation to the Maximum Influence Perturbation and its related quantities, the Most Influential Set and Perturbation Inducing Proportion, centers on a first-order Taylor expansion in $\vec{w} \mapsto \phi(\hat{\theta}(\vec{w}), \vec{w})$ around $\vec{w} = \vec{1}$. Let $\hat{\phi} := \phi(\hat{\theta}(\vec{1}), \vec{1})$, the quantity of interest at the original dataset. Then:
We can in turn approximate the Most Influential Set as follows. Let $\psi_{(n)}$ denote the order statistics of $\psi_n$, i.e., the $\psi_n$ sorted from most negative to most positive. Let $\ind{\cdot}$ denote the indicator function taking value $0$ when the argument is false and $1$ when true. Then
To compute $\vec{w}^*$ (analogous to the $\vec{w}^{**}$ that determines the exact Most Influential Set), we compute $\psi_n$ for each $n$. Then we choose $\vec{w}^*$ to have entries equal to zero at the $\lfloor \alpha N \rfloor$ indices $n$ where $\psi_n$ is most negative (and to have entries equal to one elsewhere). Analogous to the Perturbation Inducing Proportion, we can find the minimum data proportion $\alpha$ required to achieve a change of some size $\Delta$: i.e., such that $\phi^{\mathrm{lin}}(\vec{w}^*) - \hat{\phi} > \Delta$. In particular, we iteratively remove the most negative $\psi_n$ (and the index $n$) until the $\Delta$ change is achieved; if the number of removed points is $M$, the proportion we report is $\alpha = M/N$. Recall that finding the exact Maximum Influence Perturbation, Most Influential Set, and Perturbation-Inducing Proportion required running a data analysis more than $\binom{M}{\lfloor \alpha N \rfloor}$ times. By contrast, our approximation requires running just the single original data analysis, $N$ additional fast calculations to compute each $\psi_n$, and finally a sort on the $\psi_n$ values.
We define our approximate quantities, as detailed immediately above, as follows.
Below, we will sometimes emphasize that the AMIP is a sensitivity and refer to it as the AMIP sensitivity. We will say that an analysis is AMIP-non-robust if, for a particular $\alpha$ of interest, the AMIP is large enough to change the substantive conclusions of the analysis. Conversely, if the AMIP is not large enough, we say an analysis is AMIP-robust. And we generically use the AMIP acronym to describe our methodology even when calculating the Approximate Most Influential Set or Approximate Perturbation Inducing Proportion.
\seclabel{exact_lower_bound}
For any problem where performing estimation a second time is not prohibitively costly, we can re-run our analysis without the data points in the Approximate Most Influential Set and thereby provide a lower bound on the exact Maximum Influence Perturbation.
Formally, let $\vec{w}^{**}$ be the weight vector for the exact Most Influential Set, and let $\vec{w}^*$ be the weight vector for the Approximate Most Influential Set $\amis{\alpha}$. We run the estimation procedure an extra time to recover $\phi(\hat{\theta}(\vec{w}^{*}), \vec{w}^{*})$. Then, by definition,
Since $\phi(\hat{\theta}(\vec{w}^{*}), \vec{w}^{*}) - \hat{\phi}$ is a lower bound for $\mip{\alpha}$, we can use the Approximate Most Influential Set to conclusively demonstrate non-robustness. Of course, this lower bound holds for {\em any} weight vector and will be most useful if the Approximate Maximum Influence Perturbation is close to the exact Maximum Influence Perturbation. In \secref{accuracy} below, we establish the accuracy of the approximation for small $\alpha$ under mild regularity conditions.
To finish describing our approximation, it remains to detail how to compute $\psi_n = \fracat{\partial \phi(\hat{\theta}(\vec{w}), \vec{w})}{\partial w_n}{\vec{w}=\vec{1}}$ from (ref). We will refer to the quantity $\fracat{\partial \phi(\hat{\theta}(\vec{w}), \vec{w})}{\partial w_n}{\vec{w}}$ as the influence score of data point $n$ for $\phi$ at $\vec{w}$ since, as we discuss in \secref{influence_function} below, it is the empirical influence function evaluated at the datapoint $d_{n}$. To compute the influence score, we first apply the chain rule:
The derivatives of $\phi(\cdot,\cdot)$ can be calculated using automatic differentiation software baydin2017automatic,tensorflow:2015:whitepaper,jax:2018:github,pytorch:2019:lots. And once we have $\hat{\theta}(\vec{1})$ from running the original data analysis, we can evaluate these derivatives at $\vec{w} = \vec{1}$: e.g., $\fracat{\partial \phi(\theta, \vec{w})}{\partial \theta^T}{\hat{\theta}(\vec{1}), \vec{w}=\vec{1}}$.
The term $\fracat{\partial \hat{\theta}(\vec{w})}{\partial w_n}{\vec{w} = \vec{1}}$ requires slightly more work since $\hat{\theta}(\vec{w})$ is defined implicitly. We follow standard arguments from the statistics and mathematics literatures krantz2012implicit, hampel1974influence to show how to calculate it below.
Start by considering the more general setting where $\hat{\theta}(\vec{w})$ is the solution to the equation $\gamma(\hat{\theta}(\vec{w}), \vec{w}) = 0_{P} $. We assume $\gamma(\cdot, \vec{w})$ is continuously differentiable with full-rank Jacobian matrix; then the derivative $\fracat{\partial \hat{\theta}(\vec{w})}{\partial w_n}{\vec{w}}$ exists by the implicit function theorem krantz2012implicit. We can thus use the chain rule and solve for $\fracat{\partial \hat{\theta}(\vec{w})}{\partial w_n}{\vec{w}}$; in what follows, $0_{P \times N}$ is the $P \times N$ matrix of zeros.
where we can take the inverse by our full-rank assumption.
We apply the general setting above to our special case with $\gamma(\theta, \vec{w}) = \sum_{n=1}^N w_n G(\theta, d_{n})$ to find
which can again be computed with automatic differentiation software.
\seclabel{function_examples} We end this section with some concrete examples of quantities of interest. Recall from the start of \secref{MIP} that we are often interested in whether we can change the sign or significance of an estimator, or generate a significant result of the opposite sign. Recall that $\phi(\cdot)$ with only one argument is a function of $\theta$, and $\phi(\cdot, \cdot)$ with two arguments is a function of both $\theta$ and the weights $\vec{w}$.
To form our motivating examples, suppose for the remainder of this section we are interested in the $p$-th component of $\hat\theta$, where $\hat{\theta}_{p}$ is positive and statistically significant. That is, let $\hat\sigma_{p}$ be an estimator of the variance of the limiting distribution of $\sqrt{N}\hat{\phi}$, and let $\hat{\theta}_p - \frac{1.96}{\sqrt{N}} \hat\sigma_{p}$ be the lower end of our confidence interval. So we assume $\hat{\theta}_p > 0$ and $\hat{\theta}_p - \frac{1.96}{\sqrt{N}} \hat\sigma_{p} > 0$. Moreover, we will write $\hat\sigma_{p}(\theta, \vec{w})$ to emphasize that standard errors are typically given as functions of $\theta$ and the weights $\vec{w}$. For example, standard errors based on the observed Fisher information matrix $\frac{1}{N}\sum_{n=1}^N \vec{w}_n \fracat{\partial G(\theta, d_n)}{\partial \theta}{\hat{\theta}(\vec{w})}$ will, in general, depend on the weights both explicitly and through $\hat{\theta}(\vec{w})$.
To make $\hat{\theta}_{p}$ change sign, we can take
We use $-\theta_{p}$ instead of $\theta_{p}$ since we have defined $\phi$ as a function that we are trying to increase (cf.\ (ref) and the discussion after). Increasing $\phi(\hat{\theta})$, for $\phi$ in (ref), by an amount $\Delta = \hat{\theta}_{p}$ is equivalent to $\hat{\theta}_{p}$ changing sign from positive to negative.
To make $\hat{\theta}_{p}$ statistically non-significant, we wish to take the lower bound of the confidence interval to $0$. To that end, we can take
As in the previous case, we choose (ref) with a leading negative sign because we are trying to increase $\phi$ (cf.\ (ref)). Increasing $\phi(\hat{\theta}, \vec{w})$, for $\phi$ in (ref), by an amount $\Delta = \hat{\theta}_{p} - \frac{1.96}{\sqrt{N}} \hat\sigma_{p}$ is equivalent to $\hat{\theta}_{p}$ becoming statistically insignificant.
Similarly, to change to a significant result of the opposite sign, we can take
and $\Delta = \hat{\theta}_{p} + \frac{1.96}{\sqrt{N}} \hat\sigma_{p}$, for if the upper end of the confidence interval is negative, then the estimator must be negative and statistically significant.
In each case above, the quantity $\Delta$ represents how far we must move $\phi$ in order to reverse our conclusions. In this sense, $\Delta$ is a measure of the amount of “signal” in the original dataset. As we will discuss in \secref{why} below, the signal $\Delta$ is one of the three key quantities that determine AMIP robustness.
\seclabel{linear_regression}
Before continuing, we illustrate our method with an example. Economists often analyze causal relationships using linear regressions estimated via ordinary least squares (OLS), but a researcher rarely believes the conditional mean dependence is truly linear. Rather, researchers use linear regression since it allows transparent and straightforward estimation of an average treatment effect or local average treatment effect. Researchers often invoke the law of large numbers to justify the focus on the sample mean, and invoke the central limit theorem to justify the use of Gaussian confidence intervals when the sample is large. We now discuss an example from recent economics literature showing how, in practice, the omission of a very small number of data points can have outsize influence on regression parameters in the finite sample even when the full sample is large. We will study AMIP sensitivity for OLS further using simulation and theory in \secref{influence_function_ols} below.
Consider as an example the set of seven randomized controlled trials of expanding access to microcredit discussed by meager2019understanding. For illustrative purposes we single out the study with the largest sample size: angelucci2015microcredit. This study has approximately 16,500 households. A full treatment of all seven studies is in \secref{example_microcredit_linear, example_microcredit_hierarchical} along with tables and figures of the results discussed below.
We consider the headline results on household business profit regressed on an intercept and a binary variable indicating whether a household was allocated to the treatment group or to the control group. Let $Y_{ik}$ denote the profit measured for household $i$ in site $k$, and let $T_{ik}$ denote their treatment status. We estimate the following model via OLS with the regression formula $Y_{ik} \sim \beta_0 + \beta_1 T_{ik}$. In the notation of \secref{taylor_series}, we have $\theta = (\beta_0, \beta_1)^T$, $d_{ik} = (Y_{ik}, T_{ik})$ with $n = (i, k)$, and $G(\theta, d_{ik}) = (Y_{ik} - (\beta_0 + \beta_1 T_{ik})) (1, T_{ik})^T$.
We confirm the main findings of the study in estimating a non-significant average treatment effect (ATE) of -4.55 USD PPP per 2 weeks, with a standard error of 5.88. We are interested in whether we can change the sign of $\beta_1$ from negative to positive, so we take $\phi(\theta) = \beta_1$. We compute $\psi_n$ for each data point in the sample, which takes only a fraction of a second in R using our Zaminfluence package.
Examining $\vec{\psi}$, we find that one household has $\psi_n = 4.95$; removing that single household should flip the sign if the approximation is accurate. We manually remove the data point and re-run the regression, and indeed find that the ATE is now 0.4 with a standard error of 3.19. Moreover, by removing 15 households we can generate an ATE of 7.03 with a standard error of 2.55: a significant result of the opposite sign.
How is it possible for the absence of a single household to flip the sign of an estimate that was ostensibly based on all the information from a sample of 16,500? It may be tempting to suspect the use of sample means, which are known to be non-robust to gross errors, or to speculate that such excess sensitivity is simply symptomatic of ordinary sampling noise which is captured adequately by standard errors. In \secref{why} to follow, we show that such intuition is not correct. On the contrary, AMIP robustness is in fact fundamentally different than both standard errors and classical robustness to gross errors.
\seclabel{why} We now establish the determinants and accuracy of AMIP robustness. We begin by deriving the key quantities of AMIP robustness in the simple case of correctly specified univariate OLS regression (\secref{influence_function_ols}). For this simple case, we show with theory and simulations that AMIP robustness is not necessarily driven by misspecification, that AMIP non-robustness does not vanish asymptotically, and that AMIP robustness is distinct from standard errors. Next, we formally extend these conclusions to general Z-estimators in \secref{influence_function}. Finally, in \secref{accuracy}, we establish conditions under which the approximation is provably uniformly accurate for small $\alpha$, both in finite sample and asymptotically.
We will see that a central equation in our understanding of AMIP robustness is its decomposition into three key quantities: the signal, noise, and shape. First, the signal $\Delta$ is the size of change in our quantity of interest that would reverse our substantive conclusion (see \secref{function_examples} above). Large values of the signal $\Delta$ indicate that large changes are needed to make a different decision. Second, the noise $\hat{\sigma}_{\psi}$ is defined by
We call $\hat{\sigma}_{\psi}$ the noise because $\hat{\sigma}_{\psi}^2$ is typically a consistent estimator of the variance of the limiting distribution of $\sqrt{N} \phi(\hat{\theta})$, a fact that will follow below from the relationship between AMIP robustness, robust standard error estimators, and the influence function (see \secref{amip_decomposition}, paragraph (ref) or, more generally, \secref{influence_function_for_real}, paragraph (ref)). Third, the shape $\hat{\mathscr{T}}_\alpha$ is defined as
where $\psi_{(n)}$ refers to the $n$-th order statistic of the influence scores, and $\ind{\cdot}$ denotes the indicator function taking value $1$ when its argument is true and $0$ otherwise. The shape $\hat{\mathscr{T}}_\alpha$ depends in a complicated way on the shape of the tail of the distribution of the influence scores, but we show that $0 \le \hat{\mathscr{T}}_\alpha \le \sqrt{\alpha(1 - \alpha)}$ with probability one, and that $\hat{\mathscr{T}}_\alpha$ converges in probability to a nonzero constant under standard assumptions (see \secref{amip_decomposition}, paragraph (ref)). Given these three quantities, we will show in \secref{amip_decomposition}, paragraph (ref) that
We refer to the quantity $\Delta / \hat{\sigma}_{\psi}$ as the signal-to-noise ratio. For a given $\alpha$, (ref) suggests that it is the signal-to-noise ratio that primarily determines AMIP robustness. Additionally, this decomposition allows us to succinctly compare AMIP robustness to standard errors and gross-error robustness, as well as to analyze the large-$N$ behavior of AMIP robustness.
This section will use the following notation. Let the symbol $\xrightarrow{p}$ denote convergence in probability, and $\rightsquigarrow$ denote convergence in distribution, both as $N \rightarrow \infty$. Let $\vnorm{\cdot}_{op}$ denote the operator norm of a matrix.
\seclabel{influence_function_ols} We begin by focusing on the simple case of correctly-specified univariate linear regression, both to provide intuition and motivate the more general results that follow.
\point{Model} \sloppy Let $X=(x_1, \ldots, x_N)^T$ denote a vector of $N$ continuous mean-zero regressors, drawn IID from a distribution with finite variance $\sigma_x^2$. Let $\varepsilon=(\varepsilon_1, \ldots, \varepsilon_N)$ be a vector of IID draws from a $\mathcal{N}(0, \sigma_\varepsilon^2)$ distribution, where we will assume $\sigma_\varepsilon$ is known. For some unknown $\theta_0 \in \mathbb{R}$, let $y_n = \theta_0 x_n + \varepsilon_n$, so that the vector $Y=(y_1, \ldots, y_N)$ given $X$ is drawn from a correctly specified regression model with true coefficient $\theta_0$.
\point{Weighted estimating equation} The OLS estimator $\hat{\theta}$ is traditionally found by maximizing the (log) likelihood: $\log p(y_n \vert \theta, x_n) = -\frac{1}{2 \sigma_\varepsilon^{2}}(y_n - \theta x_n)^2 + C$, where $C$ does not depend on $\theta$. In particular, setting the derivative of the log likelihood to zero yields the estimating equation $G(\theta, d_n) = -\frac{1}{\sigma_\varepsilon^{2}} (y_n - \theta x_n) x_n = 0$. That is, $\hat{\theta}$ is a Z-estimator with this choice of $G$ (see (ref)). Typical Z-estimators do not have closed-form solutions. But in this case, the solution to the estimating equation returns the usual OLS estimate. A similar derivation returns the solution to the weighted estimating equation given in (ref): $\hat{\theta}(\vec{w}) = \left(\frac{1}{N}\sum_{n=1}^N \vec{w}_n x_n^2 \right)^{-1} \frac{1}{N}\sum_{n=1}^N \vec{w}_n y_n x_n$.
\point{Quantity of interest} Suppose we are interested in the sign of $\theta_0$. Without loss of generality, we assume $\hat{\theta} < 0$. Then our quantity of interest is $\phi(\theta) = \theta$.
\point{Signal and noise} For our quantity of interest, the signal is $\Delta = \@ifstar{\oldabs}{\oldabs*}{\hat{\theta}}$ since, if we can increase $\hat{\theta}$ by an amount $\@ifstar{\oldabs}{\oldabs*}{\hat{\theta}}$, its sign will change. To compute the noise, we compute the influence scores. Directly differentiating the explicit formula for $\hat{\theta}$ gives, as it must, the same value for $\psi_n$ as the implicit function theorem result of (ref). Letting $\hat\varepsilon_n := y_n - \hat{\theta} x_n$ and $S_X := \frac{1}{N}\sum_{n=1}^N x_n^2$, we see, either by direct differentiation or by (ref), that $\psi_n = N^{-1} S_X^{-1} x_n \hat\varepsilon_n$. For intuition about the noise $\hat{\sigma}_{\psi}$, we observe its asymptotic behavior. Standard results for OLS give:
Note that the noise includes a contribution from both the residual and regressor variance---we describe $\hat{\sigma}_{\psi}$ as the “noise” because it estimates the variability of $\sqrt{N} \hat{\theta}$, not of the residuals (see \secref{ols_what_determines}, paragraph (ref) below). Finally, we emphasize that, although we will be using asymptotics to provide intuition, by “noise” we will always mean the finite-sample quantity $\hat{\sigma}_{\psi}$, not its asymptotic limit.
\seclabel{ols_what_determines} Now that we have translated OLS into our framework, we can analyze the AMIP for OLS. To that end, we use both theory and a simulation study. We outline the simulation study before describing our main conclusions. For $N=5,000$ data points, and for a range of $\sigma_x$ and $\sigma_\varepsilon$, we drew normal regressors $x_n \sim \mathcal{N}(0, \sigma_x^2)$ and residuals $\varepsilon_n \sim \mathcal{N}(0, \sigma_\varepsilon^2)$. For $\theta_0 = 0.5$, we set $y_n = \theta_0 x_n + \varepsilon_n$. We computed the OLS estimator $\hat\theta = \sum_{n=1}^N y_n x_n / \sum_{n=1}^N x_n^2$.
\point{Signal-to-noise ratio drives AMIP robustness} From our discussion at the start of \secref{why}, we expect that the signal-to-noise ratio drives whether an analysis is AMIP-robust or not. In our simulation, $N$ is large and we keep $\theta_0$ fixed, so we expect that the signal does not change substantially over the simulation. Therefore, signal-to-noise is controlled by the noise. Following the asymptotic argument above, we approximate the noise as $\sigma_\varepsilon / \sigma_x$. In the left panel of \figref{sim-comb-normal}, we vary $\sigma_\varepsilon$ and $\sigma_x$ and plot the resulting Approximate Perturbation Inducing Proportion $\alpha^*$ to change the sign of $\hat\theta$. As expected, we see that the simulations with the largest approximate noise $\sigma_\varepsilon / \sigma_x$ are the least robust, in the sense that one can reverse the sign of $\hat\theta$ by dropping a very small proportion of points.
\point{Influential data points have both a large residual and large regressor} Let $(\hat\varepsilon x)_{(n)}$ denote the products $\hat\varepsilon_n x_n$, sorted from most negative to most positive, so that the sorted influence scores are $\psi_{(n)} = N^{-1} S_X^{-1} (\hat\varepsilon x)_{(n)}$. From this formula, we observe that influential datapoints have both a large residual and a large regressor (relative to the regressor variance).\footnote{ Indeed, if we had taken $\phi(\theta) = \hat{\theta} x_n = \hat{y}_n$, then the $n$-th influence score would have been $S_X^{-1} x_n^2 \hat\varepsilon_n$, which is precisely the leverage score times the residual. This expression formalizes the conceptual link made by chatterjee1986influential between influence, leverage, and large values of $\hat\varepsilon_n$. } A typical influence score goes to zero at rate $N^{-1}$, though extreme values such as $\max_{n} \@ifstar{\oldabs}{\oldabs*}{\psi_n}$ may obey a different rate. However, since $\frac{1}{N}\sum_{n=1}^N x_n^2$ and $\frac{1}{N}\sum_{n=1}^N \varepsilon_n^2$ are finite with high probability, even $\max_{n} \@ifstar{\oldabs}{\oldabs*}{\psi_n}$ does not diverge in this case.\footnote{The finiteness follows from the inequality $\frac{1}{N} \max_{n} x_n^2 \le \frac{1}{N}\sum_{n=1}^N x_n^2 \xrightarrow{p} \sigma_x^2$, with an analogous inequality for $\varepsilon_n$. However, since we know $\varepsilon_n$ is Gaussian, we actually have a stronger result in this case: $\max_{n \in \{1,\ldots,N\}} \@ifstar{\oldabs}{\oldabs*}{\varepsilon_n}$ grows at rate $\sqrt{\log (2N)}$ rigollet:2015:highdimstats.}
\point{AMIP sensitivity does not vanish as $N \rightarrow \infty$} Standard results for OLS give that $S_X \xrightarrow{p} \sigma_x^2$ and $\hat\varepsilon_n - \varepsilon_n \xrightarrow{p} 0$. So $N \psi_n - \sigma_x^{-2} x_n \varepsilon_n \xrightarrow{p} 0$. Consequently, the empirical distribution of $N \psi_n$ converges to a non-degenerate distribution with finite variance. Let $q_\alpha$ denote the $\alpha$-th quantile of the distribution of the random variable $\sigma_x^{-2} x_1 \varepsilon_1$. Since $x_n$ and $\varepsilon_n$ are independent, and about half of the $\varepsilon_n$ will be negative, we expect about half of the influence scores to be negative. So for $\alpha \ll 1/2$, with high probability at least $\alpha N$ influence scores are negative. Then, by (ref) and Slutsky's theorem, we have
The right hand side of the preceding display is strictly positive for finite $\alpha$. So, for fixed $\alpha$, we expect that AMIP sensitivity does not vanish as $N \rightarrow \infty$.\footnote{As desired, though, the expectation does go to zero as $\alpha \rightarrow 0$ since $\expect{\@ifstar{\oldabs}{\oldabs*}{x_1 \varepsilon_1}} < \infty$.}
\point{AMIP non-robustness is not due only to misspecification} Our simulations are well specified. Yet we see from \figref{sim-comb-normal} that different cases can still be robust or non-robust under various robustness cut-offs---according to their differing signal-to-noise ratios.
Asymptotically as $N \rightarrow \infty$, even in a well-specified model, we in fact expect AMIP non-robustness at any $\alpha$ for a sufficiently small $\@ifstar{\oldabs}{\oldabs*}{\theta_0}$. The limiting value of the AMIP sensitivity does not depend on $\theta_0$. Thus, as $N \rightarrow \infty$, our quantity of interest (for changing the sign of the estimator) will be AMIP non-robust with high probability if and only if $\@ifstar{\oldabs}{\oldabs*}{\theta_0} < \expect{-\frac{x_1 \varepsilon_1}{\sigma_x^2} \ind{\frac{x_1 \varepsilon_1}{\sigma_x^2} \le q_\alpha}}$. If we are interested in the sign of $\theta_0$, and $\@ifstar{\oldabs}{\oldabs*}{\theta_0}$ is small relative to the tail means of $\sigma_X^{-2} x_1 \varepsilon_1$, then the problem will be AMIP non-robust with probability approaching one, no matter how large $N$ is---despite the fact that the model is correctly specified and there are no abnormalities in the data.
\point{Though both are scaled by the noise, standard errors are different from---and typically smaller than---AMIP sensitivity} In what may seem at first like a remarkable coincidence, the variance of the limiting distribution of $N \psi_n$ (which determines AMIP sensitivity---see (ref)) is the same as the variance of the limiting distribution of our quantity of interest $\sqrt{N}(\hat{\theta} - \theta_0)$ (which determines classical standard errors). The two distributions are not the same---the limiting distribution of $N \psi_n$ is not, in general, normal---but they have the same scale. In particular, compare the noise limit in (ref) with the following limit, which follows by standard results for OLS.
As we discuss below in \secref{amip_decomposition}, paragraph (ref) and \secref{influence_function}, paragraph (ref), this equality is no coincidence, but a general (and well-known) relationship between influence scores and the limiting distributions of quantities of interest.
For large $N$, use of standard errors will admit the hypothesis that $\theta_0$ might be $0$ whenever $\@ifstar{\oldabs}{\oldabs*}{\theta_0} < \frac{1.96}{\sqrt{N}} \frac{\sigma_\varepsilon}{\sigma_x}$. Thus, for every $\theta_0 \ne 0$, using standard errors always rejects $\theta_0 = 0$ for sufficiently large $N$. By contrast, as we saw above, using the AMIP will admit a change large enough to move $\hat{\theta}$ to $0$ whenever
Thus, we see that both the AMIP sensitivity and standard errors admit larger possible values for $\hat{\theta}$ when the limiting value $\@ifstar{\oldabs}{\oldabs*}{\theta_0} / (\sigma_\varepsilon / \sigma_x)$ of the signal-to-noise ratio is large. But AMIP sensitivity is determined by the tail mean of the standardized influence scores, and standard errors are determined by a quantity that goes to zero as $N \rightarrow \infty$. Thus AMIP sensitivity is distinct from, and typically larger than, standard errors. The tail behavior of the unit-variance random variable $\frac{x_1}{\sigma_x} \frac{\varepsilon_1}{\sigma_\varepsilon}$ is exactly the shape we introduced at the start of \secref{why}. The shape captures the scale-independent shape of the tails of the distribution of the influence scores; see \secref{amip_decomposition}, paragraph (ref) below for a detailed and general analysis.
\point{Our approximation is accurate for small $\alpha$} The expression for $\hat{\theta}(\vec{w})$ depends on two terms, $\left(\frac{1}{N}\sum_{n=1}^N \vec{w}_n x_n^2 \right)^{-1}$ and $\frac{1}{N}\sum_{n=1}^N \vec{w}_n y_n x_n$, both of which are uniformly smooth functions of $\vec{w} / N$ with high probability for sufficiently small $\vnorm{\vec{w} - \vec{1}}_2 / N$. As a consequence of smoothness, we expect a linear approximation formed at $\vec{w} = \vec{1}$ to be accurate when $\vnorm{\vec{w} - \vec{1}}_2 / N$ is small. And when $\vec{w}$ contains no more than $\lfloor \alpha N \rfloor$ zeros and the rest ones, we have that $\vnorm{\vec{w} - \vec{1}}_2 / N \le \alpha$, so we expect a linear approximation to be accurate when $\alpha$ is small. We make this intuition precise and general in \secref{accuracy} below.
We check the accuracy of the approximation empirically in \figref{sim-comb-normal}. For the right hand plot in \figref{sim-comb-normal}, we fixed $\sigma_\varepsilon = 1$ and $\sigma_x = 2$. We computed the Approximate Most Influential Set for a range of left-out proportions $\alpha$ from $0$ to $10 \%$. For each $\alpha$, we computed the linear approximation, re-ran the regression to compute the actual change, and computed the error of the linear approximation as the difference of the two. The right panel of \figref{sim-comb-normal} shows how the relative error of the approximation vanishes for small $\alpha$, and that, qualitatively, the approximation is very good for removal proportions less than $2.5\%$.
\seclabel{influence_function} \def\ichat#1{\widehat{IC}_{#1}}
We next show that the conclusions of \secref{influence_function_ols} hold not just for OLS but in considerable generality for Z-estimators applied to IID data. In the present section, we will establish more generally that AMIP sensitivity is not a product of misspecification, does not vanish as $N$ goes to infinity, and is distinct from standard errors. To that end, in \secref{amip_decomposition} we first formally decompose the AMIP into the shape and noise terms defined at the beginning of \secref{why}, and we establish that the shape is roughly constant across distributions. Then, in \secref{amip_robustness_breakdown}, we use this decomposition to revisit our OLS conclusions about AMIP sensitivity but now more broadly. Finally, in \secref{influence_function_for_real}, we connect the AMIP to the influence function, showing how AMIP robustness is different from gross error robustness.
\seclabel{amip_decomposition}
\point{The AMIP is the noise times the shape} Let $\psi_{(1)}, \ldots, \psi_{(N)}$ denote the order statistics of the influence scores. Recall that the Approximate Maximum Influence Perturbation is given by the negative of the sum of the $\lfloor \alpha N \rfloor$ largest influence scores. So we can write
The first equality follows from the definition of the AMIP $\amip{\alpha}$ (\defref{approx_metrics}). The second equality follows from (ref). The third equality follows from the definitions of noise $\hat{\sigma}_{\psi}$ and shape $\hat{\mathscr{T}}_\alpha$ at the start of \secref{why}.
\point{The noise is an estimator of the standard deviation of the limiting distribution of the quantity of interest (Z-estimator version)} In the case of Z-estimators, we can show by direct computation that $\hat{\sigma}_{\psi}^2$ is the estimator of the variance of the limiting distribution of $\sqrt{N}\phi(\hat{\theta})$ given by the delta method and the “sandwich” or “robust” covariance estimator huber1967sandwich, stefanski:2002:mestimation. To see this, observe first that $\frac{1}{N}\sum_{n=1}^N \fracat{\mathrm{d} \hat{\theta}(\vec{w})}{\mathrm{d} \vec{w}_n}{\vec{1}} \left( \fracat{\mathrm{d} \hat{\theta}(\vec{w})}{\mathrm{d} \vec{w}_n}{\vec{1}}\right)^T$, as given by (ref), is precisely the sandwich covariance estimator for the covariance of the limiting distribution of $\sqrt{N} \hat{\theta}$. In turn, the sample variance of the linear approximation given in (ref), given by $\hat{\sigma}_{\psi}^2$, is then the delta method variance estimator for $\sqrt{N}\hat{\phi}$. Note that we came to the same conclusion in the special case of OLS in \secref{ols_what_determines}, paragraph (ref) above.
It follows that we can use $\hat{\sigma}_{\psi}$ to form consistent credible intervals for $\phi$, a fact that will be useful below when comparing AMIP robustness to standard errors. Specifically, if $\hat{\sigma}_{\psi} \xrightarrow{p} \sigma_{\psi}$ and $\hat{\theta} \xrightarrow{p} \theta_{\infty}$, then
As we discuss in \secref{influence_function_for_real}, paragraph (ref) below, this relationship between asymptotic variance and the influence scores is in fact a consequence of a general relationship between influence functions and distributional limits.
\point{The shape depends primarily on $\alpha$, not on the model specification} More precisely, we next show that the shape $\hat{\mathscr{T}}_\alpha$ satisfies the following properties. (1) With probability one, $0 \le \hat{\mathscr{T}}_\alpha \le \sqrt{\alpha(1-\alpha)}$. (2) Typically, $\hat{\mathscr{T}}_\alpha$ converges in probability to a nonzero constant as $N \rightarrow \infty$. (3) $\hat{\mathscr{T}}_\alpha$ is largest when the influence scores of the left-out points are all equal. Conversely, heavy tails in the distribution of $\psi_n$ result in smaller values of $\hat{\mathscr{T}}_\alpha$. (4) Empirically, $\hat{\mathscr{T}}_\alpha$ varies relatively little among common sampling distributions.
To prove the lower bound in (1), we observe that the indicator $\ind{\psi_{(n)} < 0}$ accounts for the fact that the adversarial weight would leave out fewer points rather than drop a point with positive $\psi_{(n)}$. Because of this, $\hat{\mathscr{T}}_\alpha \ge 0$. We show the upper bound of (1) as part of the extremization argument for (3) below.
To prove (2), notice that $\hat{\mathscr{T}}_\alpha$ is a sum of $\lfloor \alpha N \rfloor$ positive terms, divided by $N$. In general, then, we expect $\hat{\mathscr{T}}_\alpha$ to converge to a nonzero constant for fixed $\alpha$ as long as the distribution of $N \psi_n$ converges marginally in distribution to a non-degenerate random variable. And indeed, by (ref), we expect such convergence from Slutsky's theorem as long as $\hat{\theta}$ and $\frac{1}{N}\sum_{n=1}^N \fracat{\partial G(\hat{\theta}, d_n)}{\partial\theta}{\hat{\theta}}$ converge in probability to constants, since $N \psi_n$ is proportional to $G(\hat{\theta}, d_n)$, which itself has a non-degenerate limiting distribution.
We next show (3), that $\hat{\mathscr{T}}_\alpha$ takes its largest possible value when all the influence scores $\psi_{(1)}, \ldots, \psi_{(\alpha N)}$ take the same negative value. To that end, take $\alpha N$ to be an integer for simplicity. By the definition of $\hat{\sigma}_{\psi}$ ((ref)), $\frac{1}{N}\sum_{n=1}^N \left( \frac{N \psi_{(n)}}{\hat{\sigma}_{\psi}} \right)^2 = 1$, and by properties of the influence function detailed below, $\sum_{n=1}^N \psi_n = 0$ (\secref{influence_function_for_real}, paragraph (ref)). So $\hat{\mathscr{T}}_\alpha$ is a tail average of scalars with zero sample mean and unit sample variance. Therefore, it is equivalent to consider scalars $z_1, \ldots, z_N$ with $\frac{1}{N}\sum_{n=1}^N z_n = 0$ and $\frac{1}{N}\sum_{n=1}^N z_n^2 = 1$ and to ask how to maximize the average $-\frac{1}{\alpha N} \sum_{n=1}^{N\alpha} z_{(n)}$.
To perform this maximization we divide datapoints into a set $D$ of dropped indices, and set $K$ of kept indices. To be precise, $D := \{n: z_{(n)} \le z_{(\alpha N)} \}$ and $K := \{1,\ldots,N\} \setminus D$. We write the sample means and variances within the sets respectively as $\mu_D := \frac{1}{\alpha N} \sum_{n \in D} z_n$ and $v_D := \frac{1}{\alpha N} \sum_{n \in D} (z_n - \mu_D)^2$, with analogous expressions for $\mu_K$ and $v_K$. In this notation, our goal is to extremize $\mu_D$, the mean in the dropped set. The constraints on the distribution can then be written as $\frac{1}{N}\sum_{n=1}^N z_n = 0 \Rightarrow \alpha \mu_D + (1- \alpha) \mu_K = 0$, and $\frac{1}{N}\sum_{n=1}^N z_n^2 = 1 \Rightarrow \alpha(v_D + \mu_D^2) + (1 - \alpha) (v_K + \mu_K^2) = 1$. Given these constraints, we extremize $\mu_D$ by setting $v_K = v_D = 0$, in which case we achieve $\mu_D = -\sqrt{(1 - \alpha) / \alpha}$. Identifying $N \psi_n / \hat{\sigma}_{\psi}$ with $z_n$, and $\hat{\mathscr{T}}_\alpha$ with $\alpha \mu_D$, we see that the worst-case value of $\hat{\mathscr{T}}_\alpha$ occurs when all the influence scores $\psi_{(1)}, \ldots, \psi_{(\alpha N)}$ take the same negative value. This observation completes our argument for (3). It also follows from this argument that $\hat{\mathscr{T}}_\alpha \le \sqrt{\alpha (1 - \alpha)}$ with probability one, a bound that is achieved in the worst-case. This observation supplies the upper bound in (1).
To establish point (4), we fix a representative $\alpha$, simulate a large number of IID draws $\tilde{z}_n$ from some common distributions, standardize to get $z_n := \frac{\tilde{z}_{n} - \bar{\tilde{z}}}{\sqrt{\frac{1}{N}\sum_{n=1}^N (\tilde{z}_n - \bar{\tilde{z}})^2}}$, and compute the shape $\hat{\mathscr{T}}_\alpha = -\frac{1}{N} \sum_{n=1}^{\lfloor \alpha N \rfloor} z_{(n)}$. We find that, across common distributions, $\hat{\mathscr{T}}_\alpha$ varies relatively little. For example, for $\alpha = 0.01$, a Normal distribution gives $\hat{\mathscr{T}}_\alpha = 0.0266$, a Cauchy distribution gives $\hat{\mathscr{T}}_\alpha = 0.0022$. As expected based on the reasoning of the previous paragraph, the heavy-tailed Cauchy distribution has a smaller shape than the Normal distribution. The worst-case distribution, for which all left-out $z_n$ are equal, gives $\hat{\mathscr{T}}_\alpha = 0.0995 \approx \sqrt{\alpha(1-\alpha)}$ as expected.
\seclabel{amip_robustness_breakdown} We now use the decomposition of the AMIP into noise and shape, and the relative stability of the shape, to derive a number of general properties of AMIP robustness.
\point{Signal-to-noise ratio drives AMIP robustness} We argued above that we do not expect $\hat{\mathscr{T}}_\alpha$ to vary radically from one problem to another. By contrast, the noise $\hat{\sigma}_{\psi}$ can, in principle, be any positive number. We conclude then, that the signal-to-noise ratio, rather than the shape, principally determines AMIP robustness.
This relationship also suggests what might be done if the analysis is deemed AMIP non-robust. Since, as we showed in \secref{amip_decomposition}, paragraph (ref), $\hat{\sigma}_{\psi}$ is thus the same quantity that enters standard error computations, analysts are typically attentive to choosing estimators with $\hat{\sigma}_{\psi}$ as small as possible while still guaranteeing desirable properties like consistency. Meanwhile, the signal $\Delta$ is determined by the question being asked and the true state of nature as estimated by $\hat{\theta}$. In light of these observations, consider a case where $\Delta / \hat{\sigma}_{\psi}$ is too small to ensure AMIP robustness. Then it seems necessary for the investigator to ask a different question, or investigate different data, to find an AMIP robust analysis.
\point{AMIP sensitivity does not vanish as $N \rightarrow \infty$} Both $\hat{\sigma}_{\psi}$ and $\hat{\mathscr{T}}_\alpha$ converge to nonzero constants. So $\hat{\sigma}_{\psi} \hat{\mathscr{T}}_\alpha$, the estimated amount by which you can change an estimator, does not go to zero, either. If the signal $\Delta$ is less than the probability limit of $\hat{\sigma}_{\psi} \hat{\mathscr{T}}_\alpha$, then the problem will be AMIP non-robust no matter how large $N$ grows. As we discuss below, this behavior contrasts sharply with the behavior of standard errors.
\point{AMIP non-robustness is not due only to misspecification} Consider a correctly-specified problem with no aberrant data points. As we discussed above in \secref{amip_decomposition}, paragraph (ref), the noise will still have some non-zero probability limit. We showed in \secref{amip_decomposition}, paragraph (ref) that the shape will have a non-zero probability limit. And the quantity of interest $\phi(\hat{\theta})$ can generally be expected to have a non-zero probability limit. So by the decomposition of (ref), if the user is interested in a question whose signal is small enough, their problem will be AMIP non-robust, despite correct specification.
\point{Though both are scaled by noise, standard errors are different from---and typically smaller than---AMIP sensitivity} Recall that classical standard errors based on limiting normal approximations also depend on $\hat{\sigma}_{\psi}$, in that we typically report a confidence interval for $\phi$ of the form $\phi \in \left(\phi(\theta, \vec{1}) \pm q_{\mathcal{N}} \frac{\hat{\sigma}_{\psi}}{\sqrt{N}} \right)$, where $q_{\mathcal{N}}$ is some quantile of the normal distribution, e.g. the 0.975-th quantile $q_{\mathcal{N}} \approx 1.96$. In this sense, using standard errors errors allow that $\phi$ may be as large as $\phi + \Delta$ whenever $\Delta / \hat{\sigma}_{\psi} \le \frac{1.96}{\sqrt{N}}$. By contrast, AMIP robustness allows that $\phi$ may be as large as $\phi + \Delta$ when $\Delta / \hat{\sigma}_{\psi} \le \hat{\mathscr{T}}_\alpha$. Since $\hat{\mathscr{T}}_\alpha \ne \frac{1.96}{\sqrt{N}}$ in general, these two approaches will yield different conclusions. Indeed, typically $\hat{\mathscr{T}}_\alpha$ converges to a non-zero constant as $N \rightarrow 0$, while $\frac{1.96}{\sqrt{N}}$ converges to zero.
\point{Statistical non-significance is always AMIP-non-robust as $N \rightarrow \infty$} This observation follows as a corollary of the discussion above. In particular, we might conclude statistical non-significance if $\@ifstar{\oldabs}{\oldabs*}{\phi(\hat{\theta}, \vec{1})} \le \frac{1.96 \hat{\sigma}_{\psi}}{\sqrt{N}}$. To produce a statistically significant result, and so undermine the conclusion, it suffices to move $\phi(\hat{\theta}, \vec{1})$ by more than $\frac{1.96 \hat{\sigma}_{\psi}}{\sqrt{N}}$. Take any $\alpha$. As we have seen above, we can produce a change of $\hat{\sigma}_{\psi} \hat{\mathscr{T}}_\alpha$, which is greater than $\frac{1.96 \hat{\sigma}_{\psi}}{\sqrt{N}}$ whenever $\hat{\mathscr{T}}_\alpha > 1.96 / \sqrt{N}$. Thus, for any fixed $\alpha$, there always exists a sufficiently large $N$ such that statistical non-significance can be undermined by dropping at most $\alpha$ proportion of the data. By contrast, statistical significance can be robust if $\phi(\hat{\theta}, \vec{1})$ converges to a value sufficiently far from $0$.
\seclabel{influence_function_for_real}
We next review the influence function, its known properties, and its particular form for Z-estimators hampel1986robustbook. We first show the relationship between the influence scores and the empirical influence function. We use these connections to further justify the relationship between the noise and the limiting distribution of $\sqrt{N}\hat{\phi}$. Finally, we use these classical properties of the influence function to contrast AMIP robustness with gross error robustness and establish that outliers primarily affect AMIP robustness via the noise, rather than via the shape.
\point{Writing a statistic as a functional of the empirical distribution} Before defining the influence function, we set up some useful notation. Suppose we observe IID data, $d_1, \ldots, d_N$. Each point is drawn from a data distribution $F_{\infty}(\cdot) = p(d_1 \le \cdot)$, where the inequality may be multi-dimensional. For a generic distribution $F$, let $T$ represent a functional of the distribution: $T(F)$. One example is the sample mean; for a generic distribution $F$, let $T_{mean}(F) = \int \tilde{d} \mathrm{d} F(\tilde{d})$. Then $T_{mean}(F_{\infty}) = \expect{d_1}$ is the population mean. If we let $\hat{F}_N$ denote the empirical distribution function $\hat{F}_N(\cdot) = \frac{1}{N}\sum_{n=1}^N \ind{\cdot \le d_n}$, then $T_{mean}(\hat{F}_N) = \frac{1}{N}\sum_{n=1}^N d_n$ is the sample mean.
Now consider Z-estimators. Define $T_Z(F)$ to be a quantity satisfying
See, e.g., hampel1986robustbook. If we plug in $\hat{F}_N$ for $F$ in (ref) (and multiply both sides by $N$), we recover the Z-estimator estimating equation from (ref), with solution $\hat{\theta} = T_Z(\hat{F}_N)$. Similarly, let $\hat{F}_w$ to be the distribution function putting weight $N^{-1} w_n$ at data point $d_{n}$. Plugging in $\hat{F}_w$ for $F$ in (ref) yields the estimating equation in (ref), for weighted Z-estimators, with solution $\hat{\theta}(\vec{w}) = T_Z(\hat{F}_w)$. Finally, we can define a new functional $T_\phi(F)$ by applying the smooth function $\phi$, which picks out our quantity of interest, to $T_Z(F)$: $T_\phi(F) = \phi(T_Z(F), \vec{1})$.\footnote{ As in ordinary calculus in Euclidean space, we can also allow for explicit $F$ dependence in $\phi$ by writing $\phi(\theta, F)$. Allowing this level of generality, though, is notationally burdensome and not typical in the analysis of the influence functions for Z-estimators. So we omit this dependence for simplicity.}
\point{The influence function} The influence function $\mathrm{IF}(d; T, F)$ measures the effect on a statistic $T$ of adding an infinitesimal amount of mass at point $d$ to some base or reference data distribution $F$ reeds1976thesis,hampel1986robustbook. Let $\delta_d$ be the probability measure with an atom of size $1$ at $d$. Then
The influence function is defined in terms of an ordinary univariate derivative, and can be computed (as a function of $d$ and $F$) using standard univariate calculus. In particular, our quantity of interest has the following influence function:
By comparing (ref) with the definition of $\psi_n$ in (ref), we can see that, formally,\footnote{The factor of $N$ arises to re-write the expectation as a sum over unit-valued weights.}
\Eqref{infl_is_infl} is not a coincidence. To see this, note that the set of distributions that can be expressed as weighted empirical distributions ($\hat{F}_w$ above) is precisely the subspace of possible distribution functions concentrated on the observed data. So the derivative $N \psi_n = N \partial \phi(\hat{\theta}(\vec{w}), \vec{1}) / \partial \vec{w}_n$ ((ref)) is simply a path derivative representation of the functional derivative $\mathrm{IF}(d_n; T_\phi, \hat{F}_N)$.
We refer to the influence function applied with $F = \hat{F}_N$ as the empirical influence function hampel1986robustbook. We conclude that the $\psi_n$ that we use to form our approximation are the values of the empirical influence function at the datapoints $d_1, \ldots, d_N$. For this reason, we refer to the $\psi_n$ as influence scores.
\point{The sum of the influence scores is zero} We can now use standard properties of the influence function to reason about $\vec{\psi}$. For instance, the fact that $\sum_{n=1}^N \psi_n = 0$ follows from (ref) and the fact that $\hat{\theta}$ solves (ref).
\point{The noise is an estimator of the standard deviation of the limiting distribution of the quantity of interest (influence function version)} Observe that, by our influence function development above, we can write the squared noise as follows.
Recall that we saw above that $\hat{\sigma}_{\psi}^2$ consistently estimates the variance of the limiting distribution of $\sqrt{N}\hat{\phi}$, first in the special case of OLS (\secref{ols_what_determines}, paragraph (ref)) and then for Z-estimators in general (\secref{amip_decomposition}, paragraph (ref)). We can now see that those results are themselves special cases of the following well-known relationship between the influence function and the limiting variance of its corresponding functional:
where the expectation in the preceding display is taken with respect to $d_1 \sim F_{\infty}$ (see, e.g., hampel1986robustbook).\footnote{Though (ref) can provide useful intuition, as it does in our case, it is often easier in any particular problem to prove asymptotic results directly rather than through the functional analysis perspective of this section, since stating precise and general conditions under which (ref) holds can be challenging. See, for example, the discussion in serfling2009approximation or vaart2000asymptotic.} Specifically, if we can show that $\sigma_{\psi}$, the probability limit of $\hat{\sigma}_{\psi}$, is equal to $\expect{\mathrm{IF}(d_1; T, F_{\infty})^2}$, then (ref) would imply $\sqrt{N}(T_\phi(\hat{F}_N) - T_\phi(F_{\infty})) \rightsquigarrow \mathcal{N}(0, \sigma_{\psi}^2)$, just as we showed in (ref) using the sandwich covariance estimator. In our case, under standard assumptions, one can show directly from (ref) that $\mathrm{IF}(d_n; T_\phi, \hat{F}_N) \xrightarrow{p} \mathrm{IF}(d_n; T_\phi, F_{\infty})$, almost surely in $d_n$. A law of large numbers can then be applied to (ref) giving the desired result.
\point{AMIP robustness is different from gross error robustness} Roughly speaking, an estimator is considered non-robust to gross errors if its influence function is unbounded huber1981robust. For instance, the influence function arising from the OLS Z-estimator (\secref{influence_function_ols}) is classically known to be non-robust to gross errors. When an influence function is unbounded, one can produce arbitrarily large changes in the quantity of interest by making arbitrarily large changes to a single datapoint. Gross-error robustness is motivated by the possibility that some small number of datapoints come from a distribution arbitrarily different from the model's posited distribution. By contrast, to assess AMIP robustness, we do not make arbitrarily large changes to datapoints. We simply remove datapoints. And the analysis is AMIP-non-robust if a change of a particular size ($\Delta$) can be induced, rather than an arbitrarily large change. Consequently, problems with unbounded influence functions (such as OLS in \secref{influence_function_ols}) can be AMIP-robust if $\Delta / \hat{\sigma}_{\psi}$ is sufficiently large. And perfectly specified problems with no outliers can be AMIP non-robust if $\Delta / \hat{\sigma}_{\psi}$ is sufficiently small.
\point{Outliers affect AMIP robustness through the noise} Consideration of gross-error robustness encourages users to examine their data for unusual “outliers” in the data; once outliers are removed or their influence diminished, the problem is considered gross-error robust. Since outliers are heuristically associated with heavy-tailed data distributions, one might expect the effect of outliers to affect AMIP robustness through the shape variable $\hat{\mathscr{T}}_\alpha$. However, our analysis of \secref{amip_decomposition}, paragraph (ref) shows that gross errors actually {\em reduce} $\hat{\mathscr{T}}_\alpha$ and so render an estimator more robust for a fixed $\hat{\sigma}_{\psi}$. This observation does not imply that gross errors decrease AMIP sensitivity. Rather, gross errors increase AMIP sensitivity through the noise $\hat{\sigma}_{\psi}$. And, as we have seen, effects on $\hat{\sigma}_{\psi}$ also affect the computation of standard errors.
\seclabel{accuracy} In \secref{ols_what_determines}, paragraph (ref) we argued that our approximation was accurate in OLS for small $\alpha$. Now we extend that argument to the general case. In particular, we state sufficient conditions under which $\phi^{\mathrm{lin}}(\vec{w})$ provides a good approximation to $\phi(\hat{\theta}(\vec{w}), \vec{w})$ for small $\alpha$ uniformly for $\vec{w} \in W_\alpha$. Our key result, \thmref{theta fun_accuracy}, holds exactly in finite samples with bounds that are, in principal, computable. Additionally, the corresponding bounds can also be expected to hold with probability approaching one as $N \rightarrow \infty$ under standard assumptions.
The linear approximation we use in (ref) is a Taylor series, so its accuracy can be controlled by controlling the Taylor series residual. giordano:2019:swiss states conditions under which the first-order Taylor series approximation to $\hat{\theta}(\vec{w})$ is accurate---precisely when using the derivative as given in (ref). Under additional smoothness assumptions on $\phi$, we can extend those results to our present (ref). Since the Taylor series expansion is expressed in terms of observable non-asymptotic quantities, the resulting error bounds hold exactly in finite sample and are, in principle, computable.
We first state assumptions under which the linear approximation is accurate for the vector $\hat{\theta}(\vec{w})$.
Roughly speaking, \assuref{ij_assu} states that the estimating equation is smooth and non-singular, that the sample averages are uniformly bounded, and that the estimating equation's derivatives are Lipschitz. Other than the size of the domain $\Omega_\theta$, \assuref{ij_assu} does not depend on $W_\alpha$, nor on any asymptotic quantities; it states only (reasonable) assumptions on the actual problem at hand.
Under \assuref{ij_assu}, we are able to apply Theorem 1 of giordano:2019:swiss for $W_\alpha$ and thereby prove the uniform accuracy of a linear approximation to $\hat{\theta}(\vec{w})$ for all $\vec{w} \in W_\alpha$. To extend the accuracy of an approximation of $\hat{\theta}(\vec{w})$ to our quantity of interest $\phi$ naturally requires smoothness assumptions on $\phi$, which we now state.
We can now state our main accuracy theorem.
When $\alpha$ is small, we expect $\alpha \ll \sqrt{\alpha}$ (for example, when $\alpha = 0.01$, $\sqrt{\alpha} = 0.1 \gg 0.01$), so \thmref{thetafun_accuracy} states that the bound in the error of our linear approximation shrinks faster than the bound in the function itself as $\alpha \rightarrow 0$. In \appref{tight_bound}, we show that the dependence on $\alpha$ given in \thmref{thetafun_accuracy} are tight, i.e., that there exist problems for which the effect size and error scale as $\sqrt{\alpha}$ and $\alpha$, respectively, as $\alpha \rightarrow 0$.
\Thmref{thetafun_accuracy} is a finite-sample result, applying exactly to the problem at hand. All else equal, finite-sample results are preferable to asymptotic ones. Nevertheless, due to the many loose bounds employed in the proof, we do not expect the constants to be useful in practice. Additionally, Theorem 1 of giordano:2019:swiss may in theory require $\alpha$ to be smaller than $1 / N$, resulting in a vacuous statement. Improving these shortcomings is an important avenue for future work (e.g. giordano2019higherorder, wilson:2020:approximate). But it is therefore useful to observe that, when uniform laws of large numbers apply to $\theta \mapsto \vnorm{G(\theta, \cdot)}_2$ and $\theta \mapsto \vnorm{H(\theta, \cdot)}_2$, and the limiting functions are also non-singular, bounded, and Lipschitz, then one can expect \assuref{ij_assu} to hold with high probability and finite constants as $N \rightarrow \infty$. A precise statement of the necessary conditions for such asymptotics to apply is given in Lemma 1 of giordano:2019:swiss.
In every case we examine in our applications in \secref{examples}, we manually re-run the analysis without the data points in the removal set $\amis{\alpha}$; in doing so, we find that the change suggested by the approximation is nearly always achieved in practice (a notable exception is given and discussed at the end of \secref{example_microcredit_hierarchical}). However, linear approximations are only approximations, and intuition about the potential weaknesses of linear approximations in general apply to our approximation. The crux of \thmref{thetafun_accuracy} is that small $\alpha$ implies that $\vec{w} - \vec{1}$ is small, thus we can control the error of a linear approximation in $\vec{w}$ evaluated at $\vec{1}$. Conversely, one would not expect the approximation to work well in general for large $\alpha$ and the correspondingly larger $\vec{w} - \vec{1}$.
As an extreme example, consider when the linear approximation reports that there is no feasible way to effect a particular change; i.e., when $\aloprop{\Delta} = \texttt{NA}$ (see \defref{approx_metrics}). Such a result may seem to imply that, no matter how many datapoints one removes, the estimator will not change by an amount $\Delta$, which is often absurd. However, such a result should be taken to mean that one would have to remove such a large proportion $\alpha$ of datapoints that the linear approximation on which we are basing the $\aloprop{\Delta}$ is invalid. A more accurate interpretation of $\aloprop{\Delta} = \texttt{NA}$ is that no {\em small} proportion of points can be removed to produce a change $\Delta$, for if there were such a small proportion, the linear approximation would have discovered it.
Similarly, linear approximations cannot be expected to work well near the boundary of parameter spaces. For example, if the quantity of interest is a variance, then the true parameter is constrained to be positive, but our linear approximation is not. It can help to linearize the problem using unconstrained reparameterizations (e.g., linearly approximating the log variance rather than variance). However, as we show in \secref{example_microcredit_hierarchical}, simply transforming to an unconstrained space is still not guaranteed to produce accurate approximations near the boundary in the original, constrained space.
\seclabel{related_work} The present work belongs to an extensive “local robustness” literature, which is concerned with measuring robustness using local properties of an estimator such as series approximations. In particular, our reliance on the influence function and its related properties is shared with a great deal of the existing statistical robustness literature. Arguably beginning with mises1947asymptotic, the idea of forming series expansions in the space of data distributions was developed both for the purposes of asymptotic theory jaeckel:1972:infinitesimal,reeds1976thesis,fernholz:1983:mises,vaart:1996:empiricalprocesses, design of robust estimators hampel1974influence,hampel1986robustbook, and the detection of “outliers” (e.g., belsley:1980:regression; cain:1984:approximatecaseinfluence; cook:1986:assessment). Further, the influence function itself is a specific instance of a much broader idea of differentiating a model with respect to its inputs in order to assess sensitivity to generic perturbations (e.g. cook:1986:assessment again; diaconis:1986:bayesconsistency; ruggeri:1993:infinitesimalposteriorsensitivity; basu:1996:local, gustafson:2012:localrobustnessbook; giordano2022bnp). Our work follows in and is deeply indebted to this line of work.
The general form for the influence function of Z-estimators which we reproduce in \secref{taylor_series} has been noted many times before in the statistics literature (e.g., hampel1986robustbook; taylor:1993:unifiedapproachtoinfluentialdata; vaart2000asymptotic), the machine learning literature (e.g. koh:2017:blackbox; giordano:2019:swiss), and is of course simply a consequence of the well-known implicit function theorem krantz2012implicit. Despite this recognition, there are many examples of special cases being derived in detail for particular models (e.g. pregibon:1981:logistic; thomas:1989:assessing, hattori:2009:ipcqcasedeletion; shi:2016:gmmcasedeletion), suggesting that the simplicity of the general form of the derivative may be under-appreciated. As we argue in \secref{AMIP}, this general form is particularly useful to recognize in the age of high-quality automatic differentiation software.
Our focus on dropping data rather than “gross errors,” though not without precedent, is distinct from much of the robustness literature. Beginning with huber:1964:robustlocation, much of the statistical robustness literature has been concerned with the possibility that the model distribution may have been contaminated with an arbitrarily adversarial distribution or, equivalently, the observed dataset contains values that can take on arbitrarily misleading values. In contrast, we focus on dropping asymptotically non-vanishing amounts of data, which remains a model-agnostic data perturbation while being less adversarial --- and arguably more reasonable in certain settings, such as generalization to slightly different populations --- than data that takes on arbitrarily adversarial values.
Our “perturbation-inducing proportion” can be thought of as an example of a “breakdown point,” when the latter is defined broadly as “the proportion of data which can be changed in some way before something bad happens to the estimator.” In the tradition of concern with gross errors, the breakdown point literature is primarily concerned with the amount of data that can be changed to an arbitrary degree before an estimator can be changed by an arbitrarily large amount huber1981robust. Our concern, of course, is different: we only drop data and consider “something bad” to be a meaningful but finite change to a key quantity of interest. Early work such as huber:1983:notion raises the possibility of more generic notions of breakdown points such as ours. However, as far as the authors are aware, the present work is the first to pursue our particular notion of breakdown point in detail.
The concern with gross errors has also led to a large literature which aims to detect and define “outliers” in a context-agnostic way belsley:1980:regression,cook:1982:residualsandinfluenceinregression,cook:1986:assessment,kempthorne:1986:decision,carlin:1991:expectedutilityinfluence. Following cook:1977:detectionofinfluential, much of this literature focuses, like us, on the effect of removing datapoints, though typically only on one or a small number of datapoints. Furthermore, this line of work evaluates the effect of dropping datapoints in service of defining a context-agnostic notion of “outlier” rather than focusing, as we do, on a particular decision using the dataset at hand.
A number of authors in the outlier detection literature consider the removal of multiple points. Since their focus is always on identifying a small number of outliers, they do not consider, as we do, the inferential implications of or the accuracy of the linear approximation for leaving out a small, fixed proportion of the data. For example, hadi:2009:sensitivityinregression derives straightforward versions of classical “outlier” metrics such as Cook's distance and Andrews-Pregibon statistics for multiple datapoints. belsley:1980:regression discuss “multiple-row effects” for linear regression: motivated by the possibility that groups of points may be influential collectively but not individually, they propose a stepwise scheme for finding influential groups of observations based on repeatedly re-fitting the model, leaving the single most influential point out at each step. johnson:1983:predictiveinfluence considers the effect on a posterior predictive distribution of the removal of three points out of a set of twenty-four, which is tractable because of the closed-form solution and relatively small number of combinations. huh:1990:local observes briefly that the first-order approximation to leaving out multiple points is the sum of their influence scores, a fact which they use to produce low-dimensional visual summaries of effects of groups of observations. taylor:1993:unifiedapproachtoinfluentialdata observes that influence functions can estimate the effect of leaving out large numbers of datapoints but consider it not useful, since their primary objective is detecting small numbers of gross errors.
To the best of the authors' knowledge, our analysis of the effects of leaving out a non-vanishing proportion of the data, both on the accuracy of the empirical influence function and on inferential conclusions, is new.
\seclabel{examples}
\seclabel{example_medicaid} In our first experiment, we show that even empirical analyses that display little classical uncertainty can be sensitive to the removal of less than 1% of the sample. We consider the Oregon Medicaid study finkelstein2012oregon and focus on health outcomes. The standard errors of the treatment effects are small relative to effect size; against a null hypothesis of no effect, most $p$ values are well below 0.01. Yet we find that for most of the results, removing less than 1% of the sample can produce a significant result of the opposite sign to the full-sample analysis. In one case, removing less than 0.05% of the sample can change the significance of the result.
First we provide some context for the analysis and results of finkelstein2012oregon. In early 2008, the state of Oregon opened a waiting list for new enrollments in its Medicaid program for low-income adults. Oregon officials then drew names by lottery from the 90,000 people who signed up, and those who won the lottery could sign up for Medicaid along with any of their household members. This setup created a randomization into treatment and control groups at the household level. The finkelstein2012oregon study measures outcomes one year after the treatment group received Medicaid. About 25% of the treatment group did indeed have Medicaid coverage by the end of the trial. The main analysis investigates treatment assignment as treatment itself (“intent to treat” or ITT analysis) and uses treatment assignment as an instrumental variable for take-up of insurance coverage (“local average treatment effect” or LATE analysis).
We focus on the health outcomes of winning the Medicaid lottery, which appear in Panel B from Table 9 of finkelstein2012oregon. Each of these $J$ outcomes is denoted by $y_{ihj}$ for individual $i$ in household $h$ for outcome type $j$. The data sample to which we have access consists of survey responders ($N = 23{,}741$); some responders are from the same household. The variable $\texttt{LOTTERY}_h$ equals one if household $h$ won the Medicaid lottery, and zero otherwise. All regressions use a set of covariates $X_{ih}$ comprised of household size fixed effects, survey wave fixed effects, and the interaction between the two. All regressions also use a set of demographic and economic covariates $V_{ih}$. To infer the ITT effects of winning the Medicaid lottery, the authors estimate the following model via OLS:
To infer the LATE of taking up Medicaid on compliers, the authors employ an Instrumental Variables (IV) strategy using the lottery as an instrument for having Medicaid insurance. All standard errors are clustered on the household, and all regressions are weighted using survey weights defined by the variable weight_12m. We have access to the following seven outcome variables, presented in Panel B of Table 9 of the original paper (as well as our tables below) in the following order: a binary indicator of a self-reported measure of health being good or very good or excellent (not fair or poor), a binary indicator of self-reported health not being poor, a binary indicator of health being about the same or improving over the last six months, the number of days of good physical health in the past 30 days, the number of days on which poor physical or mental health did not impair usual activities, the number of days mental health was good in the past 30 days, and an indicator of not being depressed in last two weeks. We replicate Panel B of Table 9 of finkelstein2012oregon exactly, both for the ITT effect ($\hat{\beta}_1$) for the entire population and for the LATE on compliers ($\hat{\pi}_1$). Both analyses show strong evidence for positive effects on all health measures, with most $p$ values well below 0.01.
For each health outcome in Panel B from Table 9 of finkelstein2012oregon, we compute the AMIP to assess how many data points one needs to remove to change the sign of the treatment effect, the significance of the treatment effect, or produce a significant result of the opposite sign. The sensitivity of the LATE analysis is shown in \tableref{ohie_profit_results_iv} and the sensitivity of the ITT analysis is shown in \tableref{ohie_profit_results_reg}. In both cases we use exactly the models from the original paper, with all fixed effects and controls included and with clustering at the household level. For most outcomes, for both the LATE and ITT analysis, the sign of the treatment effect can be changed by removing around 0.5% of the data, or approximately 100 data points in a sample of approximately 22,000. The most robust outcome, “Health being better than fair” (“Health genflip 12m”), requires the removal of a little over 1% of the sample to change the sign. Across the various outcomes, we can drop even less of the sample to change the results from significant to non-significant. In some cases, we need remove only 10 or 20 data points to effect a change in significance. Finally, for most outcomes, we can remove less than 1% of the data to produce a significant result of the opposite sign. The only two exceptions, “Health genflip 12m” and “Health change flip 12m”, require the removal of slightly more than 1% to generate a significant result with the opposite sign.
We check the performance of the approximation for each analysis by re-running the model after manually removing the data points in the Approximate Most Influential Set. The result of this procedure is shown in the “Refit Estimate” column of \tableref{ohie_profit_results_iv, ohie_profit_results_reg}. For almost every result in each table, our approximate metric reliably uncovers combinations of data points that do deliver the claimed changes. As we discuss in \secref{exact_lower_bound}, the changes recorded in the “Refit Estimate” column of \tableref{ohie_profit_results_iv, ohie_profit_results_reg} form a lower bound on the true worst-case finite-sample sensitivity.
By comparing \tableref{ohie_profit_results_iv} with \tableref{ohie_profit_results_reg}, we see that the ITT results, estimated via OLS, are not notably more AMIP-robust than the LATE results, which are estimated via IV. This may seem at first counterintuitive based on a heuristic belief that IV is in some sense a less “robust” analysis than OLS in finite sample: for example, recent authors, including young2019consistency, have suggested that the uncertainty intervals for IV may be more poorly calibrated in finite samples than the intervals for OLS. However, as we discuss in \secref{why}, the quality of being “robust” in the sense of a finite-sample estimator providing a good approximation to an asymptotic quantity is simply unrelated to AMIP robustness. Neither the size of the AMIP itself nor the accuracy of the AMIP approximation depends on asymptotic arguments (see, e.g., \secref{amip_robustness_breakdown}, paragraph (ref) and the discussion of \thmref{thetafun_accuracy}). The AMIP measures the sensitivity to data ablation of a particular procedure on a particular dataset and is indifferent to the fidelity of the chosen quantity of interest to some asymptotic limit. For this reason, a procedure such as IV may be “non-robust” in the sense of having poor coverage in finite sample (as reported by young2019consistency) and yet be AMIP-robust, or vice versa. The two notions of “robustness” are simply different.
\seclabel{example_transfers} We next show that an empirical analysis can still be AMIP-non-robust even after outliers are removed. To that end, we apply our techniques to examine the robustness of the main analysis from angelucci2009indirect, one of the flagship studies showing the impact of cash transfers on ineligible (“non-poor”) households in the same villages, also known as “spillover effects.” The authors trimmed the consumption outcome for the non-poor households due to concerns about the influence of the largest values. Yet while the analysis on the poor households is quite robust, the analysis on the non-poor households---whom the trimming protocol actually affects---is much more sensitive.
angelucci2009indirect employ a randomized controlled trial to study the impact of Progresa, a social program giving cash gifts to eligible poor households in Mexico. The randomization occurs at the village level. So one can estimate both a main effect on the poor households selected to receive Progresa and also the impact on the non-eligible “non-poor” households located in the same villages as Progresa-receiving poor households.
The main results of the paper show that there are strong positive impacts of Progresa on total household consumption measured as an index both for eligible poor households and for the non-eligible households; see Table 1 of angelucci2009indirect. The variable $\texttt{C\_ind}_{it}$ denotes total household consumption for household $i$ in time period $t$. Values of $\texttt{C\_ind}_{it}$ above 10,000 are removed; such households are, by definition, non-poor. The authors study three different time periods separately to detect any change in the impact between the short and long term. They condition on a large set of variables (a household poverty index, land size, head of household gender, age, whether the household speaks an indigenous language, and literacy; at the locality level, a poverty index, and the number of households) to help ensure a fair comparison between households in the treatment and control villages. In this case these controls are important; the effects on the “non-poor” households are significant at the 5% level when the controls are included, but they are only significant at the 10% level in a simple regression on a dummy for treatment status.
The full data for the paper is available on the website of the American Economic Review thanks to the open-data policies of the journal and the authors. We can successfully replicate the results of this analysis with the controls and without, and we proceed with the controls in our present analysis in accordance with the original authors' preferred specification. We consider the time periods indexed as $t=8,9,10$ in the dataset provided, though we note that the authors do not rely on the results at $t = 8$ as the roll-out was still ongoing. We employ $K$ control variables, where $X_{itk}$ is the $k$-th variable for household $i$ in period $t$. Then we run the following regression:
Here, $\texttt{treat}_{poor,i}$ refers to an interaction between the treatment indicator and an indicator for being a poor household; correspondingly, $\texttt{treat}_{nonpoor,i}$ is an interaction between the treatment indicator and an indicator for being a non-poor household. We are able to exactly replicate the results of Table 1 of angelucci2009indirect, which exhibits positive effects of cash transfers.
We apply our methodology to assess how many data points one need remove to change the sign, the significance, or to generate a significant result of the opposite sign to that found in the full sample. We focus on the latter two time periods, as households had received only partial transfers in the first time period, but we show all three in order to replicate Table 1 from the original paper. \tableref{cash_transfers_re_run_table} shows our results. Focusing on periods 9 and 10, we find that the inferences on the direct effects on the poor households are quite robust, but the inferences on the indirect effects are less so. For the analysis of the poor, one typically needs to remove much more than 1% of the sample to change conclusions. For the analysis of the non-poor, we can remove less than 0.1% of the data to change conclusions. In fact, we can remove only 3 data points to change the significance status for both $t = 9$ and $t = 10$.
We again check the quality of our approximation. The “Refit Estimate” column in \tableref{cash_transfers_re_run_table} shows the results of manually re-running each analysis after removing the implicated data points. In most cases the AMIP correctly identifies a combination of data points that can make the claimed changes to the conclusions of the study. Although there are a few cases where re-running the analysis fails to produce the predicted statistically significant sign change, the observed changes are still large enough to be of practical interest. Furthermore, it is likely that the removal of a few additional points would in fact produce the desired statistically significant sign reversals.
Finally, we note that these results constitute an illustration of how gross error robustness is distinct from AMIP robustness (see \secref{influence_function_for_real}, paragraph (ref)). Recall that angelucci2009indirect removed (non-poor) datapoints for which consumption was greater than 10,000. By removing outliers of the consumption variable in this way, the authors of this study made what is typically considered a conservative choice in view of classical robustness concerns about gross error sensitivity. Yet, as we have shown in \tableref{cash_transfers_re_run_table}, qualitative conclusions concerning the non-poor households remain non-robust to the removal of a small number of datapoints, which demonstrates empirically that one cannot necessarily make an analysis AMIP-robust by simply trimming outliers. Indeed, as we showed above in \secref{influence_function_ols}, even perfectly specified OLS regressions with no aberrant data points can be AMIP-non-robust if the signal to noise ratio is too low.
\seclabel{example_microcredit_linear} We now show that even a simple 2-parameter linear model that performs a comparison of means between the treatment and control group of a randomized trial can be highly sensitive. To that end, we consider the analysis of seven randomized controlled trials of expanding access to microcredit, first aggregated in meager2019understanding. In \secref{example_microcredit_hierarchical} below, we will consider a more complicated Bayesian hierarchical model on the same data.
Each of the seven microcredit studies was conducted in a different country, and each study selected certain communities to randomly receive greater access to microcredit. Researchers either built a branch, or combined building a branch with some active outreach, or randomly selected borrowers among those who applied. The selected studies are: angelucci2015microcredit, attanasio2015impacts, augsburg2015impacts, banerjee2015miracle, crepon2015estimating, karlan2011microcredit, and tarozzi2015impacts. Six of these studies were published in a special issue of the American Economics Journal: Applied Economics on microcredit. All seven studies together are commonly considered to represent the most solid evidence base for understanding the impact of microcredit.
We follow the original studies and meager2019understanding in analyzing the impact of access to microcredit as the treatment of interest. The studies range in their sample sizes from around 1,000 households in Mongolia attanasio2015impacts to around 16,500 households in Mexico angelucci2015microcredit. We first focus on the headline results on household business profit regressed on an intercept and a binary variable indicating whether a household was allocated to the treatment group or to the control group. For household $i$ in site $k$, let $Y_{ik}$ denote the profit measured, and let $T_{ik}$ denote the treatment status. We estimate the following model via OLS:
This regression model compares the means in the treatment and control groups and estimates the difference as $\hat{\beta}$. We follow meager2019understanding in omitting the control variables or fixed effects from the regressions in order to examine the robustness of this fundamental procedure. But in principle this omission should make no difference to the estimate $\hat{\beta}$, and indeed it does not meager2019understanding.\footnote{The omission may in principle make a difference to the inference on $\beta$ by affecting the standard errors. However, it turns out that in these studies the additional covariates make very little difference to the standard errors. We also do not cluster the standard errors at the community level for the same reason; the results are not substantially changed. Running the regression above in each of the seven studies delivers almost identical results to the preferred specification, as it should if intra-cluster correlations are weak and covariates are not strongly predictive of household profit.}
The sensitivity results for the linear regression of profit on microcredit access appear in \tableref{mc_profit_results}. In all cases, by removing less than 1% of the data points can change either the sign or the significance. In three of the studies, one can drop less than 1% of the data points to generate a result of the opposite sign that would be deemed significant at the 5% level. Mexico, the largest study, is the most sensitive: a single data point among the 16,561 households in Mexico determines the sign (as also discussed above in \secref{linear_regression}). To produce a statistically significant result of the opposite sign---that is, to turn Mexico's noisy negative result into a “strong” positive result---one need remove only 15 data points, less than 0.1% of the sample. Mongolia, the smallest study in terms of sample size, is among the most robust in terms of sign changes; it takes 2% of the sample to change the sign. Producing a significant result of the opposite sign also requires more than 1% removal in the Philippines, Bosnia, Ethiopia, and Mongolia---whereas Mexico, India, and Morocco are more sensitive. We check the performance of our approximation by manually re-running the analysis with the data removed; the “Refit Estimate” column shows that the claimed reversal is always achieved in practice for these analyses.
By comparing the results of the present section with those of \secref{example_medicaid,example_transfers}, we can confirm the conclusion of \secref{amip_robustness_breakdown}, paragraph (ref) that standard errors are, in general, distinct from AMIP sensitivity. Despite the fact that original estimates of \tableref{mc_profit_results} are statistically insignificant, some of these non-significant results are more AMIP-robust than some of the significant results in the Cash Transfers and Oregon Medicaid examples; consider the “Significant sign change” result in the Philippines study, for example.
We next demonstrate that the AMIP sensitivity observed in \tableref{mc_profit_results} cannot simply be ascribed to statistical insignificance. To do so, we consider a different outcome with smaller variability and show that it reveals a similar sensitivity to the profit outcome. The variable we now consider is household consumption spending on temptation goods such as alcohol, chocolate, and cigarettes, since the effect of microcredit on temptation spending was estimated by meager2019understanding with the greatest precision of all six considered outcome variables. \Tableref{mc_temptation_results} shows the results of applying the AMIP to the same regression given in (ref), but with temptation spending as the outcome. While somewhat more robust than the profit analyses, the difference in the approximate removal proportions in \tableref{mc_temptation_results} is not large.
Finally, one might be tempted to ascribe the AMIP-non-robust results in \tableref{mc_profit_results} to outliers resulting from the heavy tails of the household profit variable (a phenomenon well-documented by meager2020aggregating). However, as we discuss in \secref{influence_function_for_real}, paragraph (ref) above, gross error robustness is qualitatively distinct from AMIP robustness (see also the discussion of outlier trimming at the end of \secref{example_transfers}). Indeed, the more complex hierarchical model of the next section, \secref{example_microcredit_hierarchical}, was designed precisely to accommodate the heavy tail of the household profit variable, and yet---as we will show---still exhibits a high degree of AMIP-sensitivity.
\seclabel{example_microcredit_hierarchical}
In this section, we investigate a Bayesian hierarchical model, both demonstrating that even Bayesian analyses can exhibit considerable AMIP sensitivity, and showing an example of a parameter of interest for which our linear approximation performs badly. We specifically focus on a variational Bayes approximation to the tailored mixture model from meager2020aggregating. One might hope that any of the following aspects of the more complicated model might alleviate AMIP sensitivity: the use of hierarchical Bayesian evidence aggregation, the regularization from incorporation of priors, or the somewhat more realistic data-generating process captured in this specific tailored likelihood. Indeed, the approach of meager2020aggregating was specifically motivated by the desire to capture important features of the data-generating process such as heavier tails. On the contrary, we find that the average estimated effects of microcredit remain sensitive according to the AMIP, as we did in the simpler models of \secref{example_microcredit_linear}. We also find that the linear approximation that underlies the AMIP performs poorly when attempting to decrease a particular hypervariance parameter, providing a concrete example of the limitations of our methodology, particularly for parameters near the boundary of the set of their allowable values.
Following meager2020aggregating, we fit a hierarchical model (hereafter referred to as the “microcredit model”) to all the data from the seven microcredit RCTs. We model each outcome using a spike at zero and two lognormal tail distributions, one for the positive realizations of profit and one for the negative realizations. Within the model, microcredit can affect the proportion of data assigned to each of these three components as well as affecting the location and scale of the lognormal tails. There is a hierarchical shrinkage element to the model for each parameter. The hypervariances of the treatment effects are of particular interest because these capture heterogeneity in effects across studies and offer information about the transportability of results across settings.
The models in the original paper were fit via Hamiltonian Monte Carlo (HMC) with the software package Stan carpenter2017stan. It is possible to compute the Approximate Maximum Influence Perturbation for HMC, or for any Markov Chain Monte Carlo method, using the tools of Bayesian local robustness gustafson2000local, giordano2018covariances, but the sensitivity of simulation-based estimators is beyond the scope of this paper. However, there are ways to estimate Bayesian posteriors via Z-estimators, such as with Variational Bayes (VB) techniques blei2016variational.\footnote{The Laplace approximation can also be expressed as a Z-estimator.} Specifically, we fit the microcredit model using a variant of Automatic Differentiation Variational Inference (ADVI) described in giordano2018covariances (see also the original ADVI paper, kucukelbir2017advi). Since the posterior uncertainty estimates of vanilla ADVI are notoriously inaccurate, we estimated posterior uncertainty using linear response covariances, again following giordano2018covariances.\footnote{When forming the Approximate Most Influential Set, we approximated the sensitivity only of the posterior means to data removal; the linear response covariances were considered fixed. However, when we report the results of re-fitting the model, we did re-calculate the linear response covariances at the new variational optimum.} We verified that the posterior means and covariance estimates produced by our variational procedure and the corresponding estimates from running HMC with Stan were within reasonable agreement relative to the posterior standard deviation.
We first consider the effect of microcredit on the location parameter of the positive and negative tails of profit, given respectively by the parameters $\tau_{+}$ and $\tau_{-}$. Roughly speaking, $\tau_{+}$ and $\tau_{-}$ are both estimating the effect of microcredit averaged across all of the seven countries analyzed in \secref{example_microcredit_linear}. Our point estimates for $\tau_{+}$ and $\tau_{-}$ are given by their respective VB posterior means. We used the linear response covariance estimates to form a 95% posterior credible interval in place of confidence intervals, and consider a change “significant” if the posterior credible interval does not contain zero.
\tableref{mcmix_re_run_table} shows the sensitivity of inference concerning $\tau_{+}$ and $\tau_{-}$. We see that the microcredit model's estimates of the average effectiveness of microcredit remain highly sensitive to the removal of small percentages of the sample, despite being derived from a model that accounts for non-Gaussian data shape and is regularized by the priors. This sensitivity shows that Bayesian aggregation procedures do not necessarily produce AMIP-robust estimates.
We next examine the sensitivity of the hypervariances, which measure the variability of the effect of microcredit on these tails from country to country. Specifically, the parameters $\sigma_{\tau_{+}}^2$ and $\sigma_{\tau_{-}}^2$ represent the between-country variances of the effect of microcredit on positive and negative profit outcomes, respectively. The $\sigma$ parameter can be thought of as the scale parameter analogue of the corresponding location parameter $\tau$ from \tableref{mcmix_re_run_table}. The hypervariances are of particular practical interest because they quantify how variable the effect of microcredit might be; small values of the hypervariance imply that all countries respond similarly to microcredit, whereas large values imply that one should not necessarily extrapolate the efficacy of microcredit from one country to another.
In order to avoid the possibility of extrapolating to negative variances, we form a linear approximation to our variational Bayes estimates of the posterior mean of $\log \sigma$. Since $\log \sigma$ is a scale parameter measuring the variability from country to country of the effect of microcredit, its sign is not particularly meaningful, nor is it particularly interesting to ask whether its posterior credible interval contains zero. Rather, we are interested in the magnitude of $\log \sigma$. So, to investigate robustness, we use the AMIP to check the approximate maximum change achievable in either direction (increasing or decreasing the magnitude of $\log \sigma$) by removing 0.5% of the sample, about the same fraction of the data as could generate a “significant” sign change for the $\tau_{\pm}$ parameters.
The results for the hypervariances, given in \tableref{mcmix_sd_re_run_table}, represent a useful demonstration of the limitations of our linear approximation. We are able to find sets of datapoints which, when dropped, produce {\em increases} in the hypervariances, though the our linear approximation is not nearly as accurate as in the rest of our results above. When we attempted to drop points in order to decrease the hypervariances, however, the linear approximation failed utterly; the Approximate Most Influential Set designed to produce a decrease in the hypervariances instead produced a large {\em increase} upon refitting.
Given that the hypervariances are constrained to be positive, our failure to produce large decreases may not be surprising. Note that the hypervariances' posterior expectations began very small, and that decreasing them pushes the posterior of the hypervariances closer to the boundary of the admissible space. Though the log variance may in principle take arbitrarily negative values, it nevertheless appears that the model exhibits strongly non-linear dependence on the data weights for very small variances. Designing useful diagnostics for detecting and explaining such deviations from nonlinearity in complex models is an interesting avenue for future work. In the meantime, \tableref{mcmix_sd_re_run_table} shows the importance, when possible, of checking the accuracy of the AMIP predictions by refitting the model, and of exercising caution when using the AMIP approximation near the boundary of the parameter space.
\seclabel{conclusion}
For our research conclusions to safely inform economic policy decisions, we need additional tools to quantify uncertainty beyond standard errors. There are many ways of quantifying the dependence between the finite-sample realization of the data and the conclusions of statistical inference. This dependence has become synonymous with standard errors in frequentist statistics, but the notions are equivalent only under a certain paradigm that considers a hypothetical perfect random resampling exercise for the purpose of evaluating a specific parameter within a given model. This hypothetical may not capture all the data sensitivity relevant to applied social science.
Many key ideas of 20th century statistics have their origins in the context of randomized agricultural trials, where the difference in yield across multiple fields is well-modeled by independent sampling variation. Contrast this setting with trials of economic interventions to alleviate poverty, where randomly sampling individuals or communities is a challenge and interventions may be applied across very different contexts. In such cases, statistical models are often intended to provide tractable and interpretable summaries or proxies of the general impact of interventions. As such methods often average information across individuals in ways that may not always reflect broader policy interests, it seems essential to interrogate the sensitivity of our conclusions to departures from the hypothetical thought experiment.\footnote{In agricultural trials, total yields are the true quantity of interest; for microcredit trials, the average treatment effect is but a convenient summary. If the average profit were to increase slightly through one individual becoming wealthy while leaving all others destitute, one could consider the intervention a failure. By contrast, if a single plant produced an entire harvest's worth of corn, the outcome would still be desirable, if strange.}
In this paper, we have offered one alternative way of conceiving of and quantifying the dependence of empirical results on the sample data, beyond standard errors. Sensitivity of conclusions to data removal under our metric does not necessarily imply a problem with the sample. But the goal of inference is not to learn about the sample, but to learn about the population. If minor alterations to the sample can generate major changes in the inference, and we know that the environment in which we do economics is changing all the time, we ought to be less confident that we have learned something fundamental about the world we seek to understand, for which we ultimately seek to make policy. This does not imply that the original analysis is invalid according to classical sampling theory, and we do not recommend that researchers abandon the original full-sample results even if they are sensitive according to our metric. However, reporting our metrics alongside standard errors would improve our ability to understand and interpret the findings of a given analysis.
Since AMIP analysis always indicates which data points have high (approximate) influence, our methods allow researchers not only the chance to check that the approximation worked on their own sample, but to understand what---if anything---makes these data points special. Investigating influential points may provide insight into the way in which a given inferential procedure is using the finite-sample information to generate claims about the population parameters. In addition, in cases when this sensitivity is undesirable, it may be fruitful to develop new statistical methods to ameliorate it.