# What is the name for this unintuitive result with OLS on a very “asymmetric” regressand, and how should it be addressed?

Say our sample consists of about a hundred Belgian (x = 0) and Swiss (x = 1) chocolate bars. We test them to see if they have safe (y = 1) or lethal (y = 0) levels of arsenic. As it turn out, 90% are safe. We then regress y on x with OLS to check whether the country of origin is correlated with safety:

------------------------------------------------------------------------------
y |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
swiss |        -.9   .3030152    -2.97   0.004    -1.501248   -.2987522
_cons |         .9   .0301511    29.85   0.000     .8401736    .9598264
------------------------------------------------------------------------------


Goodness! I'm never going to eat a Swiss chocolate bar again. The point estimate is exceedingly large, and the test is exceedingly significant. But looking at the data more closely, there is only one Swiss bar in our sample. It so happens that it is unsafe.

In other words, if the null hypothesis is true (i.e. x is uncorrelated with y), and if we were to repeat this experiment arbitrarily many times with a single “x = 1” observation in each sample, then 10% of the time our model would lead us to reject the null with 99.6% confidence (i.e. a type I error).

1. Does this phenomenon have a name?
2. Is there an intuitive explanation for what appears to be a gross predilection to type I errors?
3. Should this result be "corrected"? If so, how? (Without collecting more data, of course.)
4. Is there a way to identify this phenomenon in published tables, if an accompanying summary table of descriptive statistics is not presented? Are there any red flags in the regression results below?

Stata code to replicate this result:

// Create data (101 obs)
clear
input   y x n
0 1 1
0 0 10
1 0 90
end
expand n

// Describe data
tab y x

// OLS regression
reg y x

/*
// A variety of alternative estimations, with SE listed
reg y x, vce(ols)       // 0.3
reg y x, vce(robust)    // 0.03
reg y x, vce(cluster x) // ~0
reg y x, vce(bootstrap) // 0.03
reg y x, vce(jackknife) // 0.03
reg y x, vce(hc2)       // 0.03
reg y x, vce(hc3)       // 1.1
*/

exit


Results:

.     // Describe data
.     tab y x

|           x
y |         0          1 |     Total
-----------+----------------------+----------
0 |        10          1 |        11
1 |        90          0 |        90
-----------+----------------------+----------
Total |       100          1 |       101

.     // OLS regression
.     reg y x

Source |       SS           df       MS      Number of obs   =       101
-------------+----------------------------------   F(1, 99)        =      8.82
Model |  .801980198         1  .801980198   Prob > F        =    0.0037
Residual |           9        99  .090909091   R-squared       =    0.0818
Total |   9.8019802       100  .098019802   Root MSE        =    .30151

------------------------------------------------------------------------------
y |      Coef.   Std. Err.      t    P>|t|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
x |        -.9   .3030152    -2.97   0.004    -1.501248   -.2987522
_cons |         .9   .0301511    29.85   0.000     .8401736    .9598264
------------------------------------------------------------------------------

• Exactly what is this "phenomenon" to which you refer? How do you deduce that "if we were to repeat this experiment arbitrarily many times with a single “x = 1” observation in each sample, then 10% of the time our model would lead us to reject the null with 99.6% confidence (i.e. a type I error)." It's not even clear what "this experiment" is. – whuber May 5 '15 at 23:29
• You sample consists of 101 items. – Alecos Papadopoulos May 6 '15 at 0:59
• Red flags: grossly unbalanced sample and also zero in the count table. – Alecos Papadopoulos May 6 '15 at 1:19
• Alecos: (1) Alright, I've changed "a hundred" to "about a hundred". (2) But how would I know if the sample were grossly unbalanced if the count table were not present? – Nils Enevoldsen May 6 '15 at 2:31
• Side note: you are modelling a binary (0,1) outcome using a linear regression model; this does not sound very good. However, this comment does not answer your main question. – Richard Hardy May 6 '15 at 6:36

@RichardHardy's observation that you're modelling a dichotomous response using ordinary least squares regression is at the heart of the matter: the distribution of errors is far from normal, & with a single observation at $x=1$, the distribution under the null of the standardized coefficient estimate for $x$ will show a gross departure from a t-distribution.

Your description of the situation & the predictor distribution are enough to raise a red flag, but the usual regression diagnostics will highlight the issue: a clearly non-normal distribution of the residuals & the hugely disproportionate leverage of a single point. In general tables of coefficient estimates, standard errors, t-statistics, &c. don't help you to check whether the model's badly mis-specified (that's why we have diagnostics)—in this case, knowing that $x$ is dichotomous & that reference-level coding was used, & without other predictors to muddy the waters, the large difference in standard errors between the estimates for the intercept & for the $x$ coefficient reflects the discrepancy in sample size between the two groups.

# make data & fit model
x <- c(rep(0,100),1)
y <- c(rep(1,90),rep(0,10),0)
lm(y~x) -> lmod

# leverage plot
plot(hatvalues(lmod), type="h")

# normal quantile-quantile plot
qqnorm(residuals(lmod))

# simulate y values under the null of common 10% probability for Y=1
t.stat <- numeric(10000)
for (i in 1: 10000){
y.sim <- rbinom(101, 1, 0.9)
lm(y.sim~x) -> mod.sim
summary(mod.sim)$coef["x", "t value"] -> t.stat[i] } # plot kernel-smoothed density estimate of t-statistic distribution # for the x coefficient plot(density(t.stat)) # plot empirical distribution function of t-statistic plot(ecdf(t.stat)) # compare to theoretical distribution t.cumve <- function(x) pt(x, df=lmod$df.residual)