# Practical limits to collinearity problems?

Collinear independent variables can have undesirable effects on the interpretation of coefficients in a linear model. Indeed, for two perfectly-correlated predictors, the coefficients are not uniquely determined, leaving a single degree of freedom by which they may vary.

However, adding in even a little bit of randomness gives solutions to the OLS equations that are not too far from assigning the correlated independent variables equal weight.

An example. The following (python) code creates a simple linear relationship between a single independent variable and a single dependent variable. It then copies the independent variable to generate a second independent variable and adds a tiny amount of independent gaussian noise to both. Thus, the second variable is nearly-perfectly correlated with the first variable. When we do this a number of times and plot the resulting coefficients, they tend to center around 2.5 and 2.5:

import numpy as np
import matplotlib.pyplot as plt
import sklearn.linear_model

coefs = []                                # to hold the coefficient of all of the OLS fits
for i in range(1000):                     # run OLS a bunch to see what the coefficients do
X = np.linspace(0,1,100)              # some independent variable
y = 5*X + 4                           # some dependent variable
X2 = X                                # a new independent variable that is correlated with the first
XX = (
np.stack((X,X2)).T +              # stack the independent variables
)
lr = sklearn.linear_model.LinearRegression()
lr.fit(XX,y)
coefs.append(lr.coef_)                # get the coefficients of an OLS linear regression

coefs = np.array(coefs)
plt.scatter(coefs[:,0], coefs[:,1])
plt.show()


Why am I relatively unlikely to get, say, an X1 coefficient of -105 and an X2 coefficient of 110? Those add up to 5, but there is something pushing the results toward 2.5, 2.5. What is the intuition behind this phenomenon, and what practical implications does this have when faced with collinear independent variables?

UPDATE:

Altering the amount of Gaussian noise to have a much smaller variance ($$\sigma = 0.00000001$$ as opposed to $$\sigma=0.01$$ in the example above):

And $$\sigma=0.0000000000000001$$:

UPDATEUPDATE:

Interest in how the amount of added noise plays into this phenomenon encouraged me to make the following graph. I ran the above experiment for various quantities of added noise, from $$2^{1}$$ to $$2^{-60}$$ on an exponential scale. I then measured the std of the X1 coefficient in the above experiment for each of those trials. Results are given here (should read $$-\log_2$$ on the x-axis):

I suspect that at very very small values of added noise, there are numerical problems due to limits of the computer architecture. Accounting for that, it seems from this graph that the distribution of OLS coefficients for X1 and X2 converges to a Gaussian distribution with mean 2.5 and std 0.25 or so. Thus, this seems like "Why not coefficients (-105,110)?" is still open.

Zooming in on the weird behavior for very small values:

Which seems like it's doing a sane thing. That doesn't mean that it's not floating-point errors, but it's not clear that that is what's happening.

• Because I don't know the details of your programming environment, I cannot determine what you mean by a "tiny" amount of noise. How much is it and on what basis do you determine it is "tiny"? You are implicitly comparing the distribution of fitted values to some other distribution (which presumably has a greater spread), but what would that other distribution be? What were you expecting? (It's still an interesting question, btw.) – whuber Sep 27 '18 at 18:40
• I suspect that "tiny" is highly relevant. The study you are beginning to carry out in your comment concerns how the range of estimated coefficients depends on the variance of the noise in the data. It would be interesting to pursue that further, because in so doing you might gain a great deal of insight into the question and possibly even develop an answer. I also suspect what you see in your experiment may depend fairly strongly on how many iterations you are performing, so it would be good to visualize the results in a more perspicuous fashion. – whuber Sep 27 '18 at 18:47
• Ah. Seems you are correct. Setting $\sigma=0.0000000000000001$ gets the coefficients out to about $(0,5)$. Still, for such a small amount of randomness to constrain the possible OLS solutions to such a degree strikes me as odd. – Scott Sep 27 '18 at 18:50
• Yes, undertaking a quantitative study of this effect makes the question more interesting, doesn't it? :-) – whuber Sep 27 '18 at 18:52
• You will learn more by (a) varying the coefficient of X in the definition of y and (b) multiplying X2 by different factors. I believe the initial rise to a constant sd is an artifact of your software: it does not occur in R's lm function, for instance. This latter fact is a strong clue about what's going on. – whuber Sep 27 '18 at 20:12

### Why am I relatively unlikely to get, say, an x_1 coefficient of -105 and an x_2 coefficient of 110? Those add up to 5, but there is something pushing the results toward 2.5, 2.5.

Linear combinations of your two variables $$x_1 = x + \epsilon_1$$ and $$x_2 = x+\epsilon_2$$ can be described like:

$$\frac{a+b}{2} (x + \epsilon_1) + \frac{a-b}{2}(x + \epsilon_2) = a x + \frac{1}{2} a (\epsilon_1 + \epsilon_2) + \frac{1}{2} b (\epsilon_1 - \epsilon_2)$$

• The parameter a will be approximately equal to the parameter associated with the variale $$x$$. In your case that is $$a = 5$$.
• The parameter b will be related to the variance and correlation of $$y=(\epsilon_1 + \epsilon_2)$$ and $$z=(\epsilon_1 - \epsilon_2)$$ by: $$\text{Var} \left( a y + b z \right) = a^2 \text{Var} (y) + b^2 \text{Var}(z) + 2ab \sqrt{\text{Var}(y)\text{Var}(z)} \rho_{y,z}$$ note that $$y=(\epsilon_1 + \epsilon_2)$$ and $$z=(\epsilon_1 - \epsilon_2)$$ are iid distribtued variables and $$\rho_{y,z}$$ will be distributed around zero.

So mostly you will be close to $$b=0$$

So in the case of -105 and 110 you would get larger contributions from the error terms, which only gets 'undone' when there is an a strong correlation in the particular sample of the error terms.

### Influence of $$\sigma$$

I can model as well the influence of $$\sigma$$ but I do not get the same pattern as you. Below you see that with larger variance, the sum of the parameters will be smaller than $$5$$ (to decrease the effect of the error terms) and also the difference will be smaller as it relates to the size of the sum of the parameters. But, I do not see why the parameters would go from 0 to 5 as in your last graph.

This is for thousand repetitions of the data:

$$\begin{array}{rcl} x_1 &=& x + \epsilon_1 \\ x_2 &=& x + \epsilon_2 \\ y &=& 5x + 4 \end{array}$$

where $$x$$ is a vector of size $$n=100$$ varying from 0 to 1, $$\epsilon_1$$ and $$\epsilon_2$$ are Gaussian random noise.

Which is modeled as a linear model minimizing the sum of least squared error $$\epsilon$$

$$y = a x_1 + b x_2 + \epsilon$$

In the way you've added noise, you could write $$x' = x +\epsilon$$ (where $$\epsilon$$ is a normally distributed variable representing noise). Furthermore, $$y=5x+4+\eta$$, where $$\eta$$ is another Gaussian noise term.

You are trying to fit a regression of the form $$Ax + Bx' +C$$, and you know that you're fitting it to a target variable which was generated by $$y=5x+4+\eta$$, thus :

$$Ax + B(x + \epsilon) + C= 5x+4+\eta$$

or

$$(A+B)x + C + B\epsilon = 5x + 4 + \eta$$

Hopefully, this sheds some light on what's happening (together with your numerics), even if it doesn't formally show it. If $$\epsilon$$ is so small that $$B\epsilon$$ is also very small, then pretty much any combination of (A,B) which satisfies $$A+B=5$$ will likely be a good fit, and which one is best is a question of noise. The large $$\epsilon$$ is, the more $$B\epsilon$$ will in general make the LHS fluctuate wrt to the RHS and the more B will go to zero ( a trivial restatement of the fact that x is the truly explanatory variable and B is only correlated, and as you shrink that correlation, the regression will find it easier to learn this)

It's less clear to me what happens when $$\epsilon$$ and $$\eta$$ are of similar sizes. In general there will be competing effects where you might be able to match the data better with non-zero B. What happens if you simulate a larger data set? Does the effect still hold?

• I'm not adding a $\eta$, and the two copies of $x$ have i.i.d gaussian noise. Maybe $(A+B)x + C + A\epsilon_1+B\epsilon_2 = 5x+4$. – Scott Sep 27 '18 at 19:49
• I think, from looking at your code, that that's not quite right. Firstly, you set $y=5x+4$, and then you add noise to x, meaning that $y = 5(x-\epsilon _{1})+4$ and $x' = x + \epsilon_{2}$. Combining this, you get: $Ax + Bx' + C = Ax + B(x + \epsilon_{2})+C= (A+B)x + B\epsilon_{2}+C=5x - 5\epsilon_{1}+4$, which is equivalent to what I wrote above, but perhaps more insightful, because you see the relative size of the error terms – gazza89 Sep 28 '18 at 9:09
• apologies, I think it's actually slightly different. You set x first, and then y=5x+4, but then you add $\epsilon_{1}$ to x, so in terms of the current value of x, $y=5(x-\epsilon_{1})+4$. Also, x' is equal to the old value of x (pre noise) + $\epsilon_{2}$, so $x'=x-\epsilon_{1}+\epsilon_{2}$. Putting all this together, $(A+B)x + B(\epsilon_{2}-\epsilon_{1})+C = 5x-5\epsilon_{1}+4$ – gazza89 Sep 28 '18 at 10:14