# When would least squares be a bad idea?

If I have a regression model: $$Y = X\beta + \varepsilon$$ where $\mathbb{V}[\varepsilon] = Id \in \mathcal{R} ^{n \times n}$ and $\mathbb{E}[\varepsilon]=(0, \ldots , 0)$,

when would using $\beta_{\text{OLS}}$, the ordinary least squares estimator of $\beta$, be a poor choice for an estimator?

I am trying to figure out an example were least squares works poorly. So I am looking for a distribution of the errors that satisfies previous hypothesis, but yields bad results. If the family of the distribution would be determined by mean and variance that would be great. If not, it's OK too.

I know that "bad results" is a little vague, but I think the idea is understandable.

Just to avoid confusions, I know least squares are not optimal, and that there are better estimators like ridge regression. But that's not what I am aiming at. I want an example were least squares would be unnatural.

I can imagine things like, the error vector $\epsilon$ lives in a non-convex region of $\mathbb{R}^n$, but I'm not sure about that.

Edit 1: As an idea to help an answer (which I can't figure how to take further). $\beta_{\text{OLS}}$ is BLUE. So it might help to think about when a linear unbiased estimator would not be a good idea.

Edit 2: As Brian pointed out, if $XX'$ is bad conditioned, then $\beta_{\text{OLS}}$ is a bad idea because variance is too big, and Ridge Regression should be used instead. I'm more interested is in knowing what distribution should $\varepsilon$ in order to make least squares work bad.

$\beta_{\text{OLS}} \sim \beta+(X'X)^{-1}X'\varepsilon$ Is there a distribution with zero mean and identity variance matrix for $\varepsilon$ that makes this estimator not efficient?

• I don't want to sound harsh, but I'm not entirely sure what you want. There are lots of ways something could be a poor choice. Typically, we evaluate estimators in terms of things like bias, variance, robustness, and efficiency. Eg, as you note, the OLS estimator is BLUE. – gung Jan 21 '14 at 1:48
• OTOH, the variance might be so large as to be useless, making a lower-variance but biased estimator like ridge preferable. Another example is that OLS maximally uses all of the info in your data, but this makes it susceptible to outliers. There are lots of alternative loss functions that are more robust, while attempting to maintain efficiency. It might be clearer if you could re-frame your question in terms like these. I don't know what it means for an estimator to be "unnatural". – gung Jan 21 '14 at 1:49
• Thanks for your comment, it made me realize the ambiguity of the question. I hope it's clearer now – Manuel Jan 21 '14 at 16:59
• See the regression in this answer. In short: influential outliers can be a problem. – Glen_b Jan 30 '15 at 11:19

Brian Borchers answer is quite good---data which contain weird outliers are often not well-analyzed by OLS. I am just going to expand on this by adding a picture, a Monte Carlo, and some R code.

Consider a very simple regression model: \begin{align} Y_i &= \beta_1 x_i + \epsilon_i\\~\\ \epsilon_i &= \left\{\begin{array}{rcl} N(0,0.04) &w.p. &0.999\\ 31 &w.p. &0.0005\\ -31 &w.p. &0.0005 \end{array} \right. \end{align}

This model conforms to your setup with a slope coefficient of 1.

The attached plot shows a dataset consisting of 100 observations on this model, with the x variable running from 0 to 1. In the plotted dataset, there is one draw on the error which comes up with an outlier value (+31 in this case). Also plotted are the OLS regression line in blue and the least absolute deviations regression line in red. Notice how OLS but not LAD is distorted by the outlier:

We can verify this by doing a Monte Carlo. In the Monte Carlo, I generate a dataset of 100 observations using the same $x$ and an $\epsilon$ with the above distribution 10,000 times. In those 10,000 replications, we will not get an outlier in the vast majority. But in a few we will get an outlier, and it will screw up OLS but not LAD each time. The R code below runs the Monte Carlo. Here are the results for the slope coefficients:

               Mean   Std Dev   Minimum   Maximum
Slope by OLS   1.00      0.34     -1.76      3.89
Slope by LAD   1.00      0.09      0.66      1.36


Both OLS and LAD produce unbiased estimators (the slopes are both 1.00 on average over the 10,000 replications). OLS produces an estimator with a much higher standard deviation, though, 0.34 vs 0.09. Thus, OLS is not best/most efficient among unbiased estimators, here. It's still BLUE, of course, but LAD is not linear, so there is no contradiction. Notice the wild errors OLS can make in the Min and Max column. Not so LAD.

Here is the R code for both the graph and the Monte Carlo:

# This program written in response to a Cross Validated question

# The program runs a monte carlo to demonstrate that, in the presence of outliers,
# OLS may be a poor estimation method, even though it is BLUE.

library(quantreg)
library(plyr)

# Make a single 100 obs linear regression dataset with unusual error distribution
# Naturally, I played around with the seed to get a dataset which has one outlier
# data point.

set.seed(34543)

# First generate the unusual error term, a mixture of three components
e <- sqrt(0.04)*rnorm(100)
mixture <- runif(100)
e[mixture>0.9995] <- 31
e[mixture<0.0005] <- -31

summary(mixture)
summary(e)

# Regression model with beta=1
x <- 1:100 / 100
y <- x + e

# ols regression run on this dataset
reg1 <- lm(y~x)
summary(reg1)

# least absolute deviations run on this dataset
reg2 <- rq(y~x)
summary(reg2)

# plot, noticing how much the outlier effects ols and how little
plot(y~x)
abline(reg1,col="blue",lwd=2)
abline(reg2,col="red",lwd=2)

# Let's do a little Monte Carlo, evaluating the estimator of the slope.
# 10,000 replications, each of a dataset with 100 observations
# To do this, I make a y vector and an x vector each one 1,000,000
# observations tall.  The replications are groups of 100 in the data frame,
# so replication 1 is elements 1,2,...,100 in the data frame and replication
# 2 is 101,102,...,200.  Etc.
set.seed(2345432)
e <- sqrt(0.04)*rnorm(1000000)
mixture <- runif(1000000)
e[mixture>0.9995] <- 31
e[mixture<0.0005] <- -31
var(e)
sum(e > 30)
sum(e < -30)
rm(mixture)

x <- rep(1:100 / 100, times=10000)
y <- x + e
replication <- trunc(0:999999 / 100) + 1
mc.df <- data.frame(y,x,replication)

ols.slopes <- ddply(mc.df,.(replication),
function(df) coef(lm(y~x,data=df))[2])
names(ols.slopes)[2] <- "estimate"

function(df) coef(rq(y~x,data=df))[2])

summary(ols.slopes)
sd(ols.slopes$estimate) summary(lad.slopes) sd(lad.slopes$estimate)

• @Manuel Thanks. I caught an error in my R program---there was a 0.04 where there should have been an sqrt(0.04). It did not change the thrust of the answer. It made a small difference to the results. However, if you copied the code before, you should copy it again now. – Bill Jan 21 '14 at 19:08

One example would be where you do not want to estimate the mean. This came up in work I used to do where we were estimating the number of sex partners people had, as part of modelling the spread of HIV/AIDS. There was more interest in the tails of the distribution: Which people have many many partners?

In this case, you could want quantile regression; an underused method, in my opinion.

• what do you mean with not wanting to estimate the mean? I am considering $\beta$ as a fixed parameter in a frequentist approach, if that what you are talking about. – Manuel Jan 20 '14 at 23:05
• Or did you meant, the mean of $Y$ ? – Manuel Jan 20 '14 at 23:08
• Yes, I meant the mean of Y. That is what OLS regression does. – Peter Flom Jan 20 '14 at 23:30
• Great point,+1. Just consider the expression $(x-a)^{2}+(x-b)^{2}$ and calculate its minimum. It is the midpoint between a and b. This is contrast with the $L_{1}$ error function, a.k.a. robust regression – jpmuc Jan 20 '14 at 23:32

If $X$ is a badly conditioned matrix or exactly singular, then your least squares estimator will be extremely unstable and useless in practice.

If you limit your attention to the distribution of $\epsilon$, then you should keep in mind that the Gauss-Markov theorem ensures that the least squares solution will be a minimum variance unbiased estimator.

However, if the distribution of $\epsilon$ is sufficiently extreme, then it's possible to construct examples where the distribution of the estimates has bad properties (in particular, the possibility (albeit with low probability) of extremely large errors in $\beta$) despite being minimum variance.

• What would be an extreme distribution of $\varepsilon$ ? Remember it has identity covariance matrix. – Manuel Jan 20 '14 at 23:19
• Besides, Gauss-Markov ensure least squares be a minimum variance unbiased estimator among linear. May be linear estimators are not reazonable for some kind of distributions. That's what i want to understand. – Manuel Jan 20 '14 at 23:20
• Consider the distribution where $\epsilon_{i}=0$ with probability 0.9999, and $\epsilon_{i}=100$ with probability $0.00005$ and $\epsilon_{i}=-100$ with probability $0.00005$. Now, suppose $X=I$ (the y's are just direct measurements of the unknown parameter $\beta$), and you have about 100 observations. It's most likely that your estimate will be perfect, but there's a significant possibility of an estimate that includes one of the rare bad $\epsilon$ values and is off as a result. – Brian Borchers Jan 20 '14 at 23:26
• If you want something really weird, consider using an appropriately scaled Student's t distribution with 4 degrees of freedom. This is a well known distribution with finite mean and variance but unbounded fourth moment. Now, suppose that you have 1 observation with $X=1$, and $\beta=0$. The distribution of $\hat{\beta}$ will have finite mean and variance but unbounded fourth moment. – Brian Borchers Jan 20 '14 at 23:29
• I belive that your second comment is what i am looking for. It just bothers me the low probability of happening. Also from that example it's pretty clear how to build a better estimator when knowing the distribution of the errors. – Manuel Jan 21 '14 at 17:02