Why does pre-summing variables lead to a different OLS fit? Suppose I have an OLS model like this:
$$y = \beta_1x_1 + \beta_2x_2 + \epsilon$$
If you sum the variables first, $x_3 = x_1 + x_2$ and fit
$$y = \beta_3x_3 +\epsilon$$
I expect that $\beta_3$ is the mean of $\beta_1$ and $\beta_2$. In fact it's very close but not exact. The estimated intercept is also slightly different. Why is that?
x1 <- rnorm(100)
x2 <- rnorm(100)
x <- x1 + x2

y <- x1 + x2 + rnorm(100)

m0 <- lm(y ~ x1 + x2)
m1 <- lm(y ~ x)

mean(coef(m0)[2:3]) - coef(m1)[2]
```

 A: You would only get the mean if $x_1$ and $x_2$ are normalized and uncorrelated.
To see this:
If we let $\vec{x_3} = \vec{x_1} + \vec{x_2}$, we get the weight $$w_3 = \Big(\vec{x_1}^T\vec{x_1} + 2\vec{x_1}^T\vec{x_2} + \vec{x_2}^T\vec{x_2}\Big)^{-1}(\vec{x_1} + \vec{x_2})^T\cdot\vec{y} = \frac{(\vec{x_1} + \vec{x_2})^T\cdot\vec{y}}{\vec{x_1}^T\vec{x_1} + 2\vec{x_1}^T\vec{x_2} + \vec{x_2}^T\vec{x_2}}$$.
If we compare that with using both features in a combined matrix $X = \Bigg(\vec{x_1}\hspace{4mm} \vec{x_2}\Bigg)$, there we would get the weight vector$$\vec{w} = \frac{1}{D} \begin{pmatrix} (\vec{x_2}^T\vec{x_2})(\vec{x_1}^T\vec{y}) - (\vec{x_1}^T\vec{x_2})(\vec{x_2}^T \vec{y}) \\ -(\vec{x_1}^T\vec{x_2})(\vec{x_1}^T\vec{y}) + (\vec{x_1}^T\vec{x_1})(\vec{x_2}^T\vec{y})\end{pmatrix}$$ where $D$ is the determinant of the covariance matrix $X^TX$, so $D = (\vec{x_1}^T\vec{x_1})(\vec{x_2}^2\vec{x_2}) - (\vec{x_1}^T\vec{x_2})^2$.
Now, let's assume $\vec{x_1}^T\vec{x_2} = 0$ (i.e. the variables are independent). This makes $$w_3 = \frac{(\vec{x_1} + \vec{x_2})^T\cdot\vec{y}}{\vec{x_1}^T\vec{x_1} + \vec{x_2}^T\vec{x_2}}$$ and $$\vec{w} = \frac{1}{(\vec{x_1}^T\vec{x_1})(\vec{x_2}^2\vec{x_2})} \begin{pmatrix} (\vec{x_2}^T\vec{x_2})(\vec{x_1}^T\vec{y}) \\ (\vec{x_1}^T\vec{x_1})(\vec{x_2}^T\vec{y})\end{pmatrix}  = \begin{pmatrix} (\vec{x_1}^T\vec{y})/(\vec{x_1}^T\vec{x_1}) \\ (\vec{x_2}^T\vec{y})/(\vec{x_2}^T\vec{x_2})\end{pmatrix} $$
The mean of the two weights would be $\frac{(\vec{x_2}^T\vec{x_2})\vec{x_1}^T + (\vec{x_1}^T\vec{x_1})\vec{x_2}^T}{2(\vec{x_1}^T\vec{x_1})(\vec{x_2}^T\vec{x_2})}\cdot\vec{y}$.
So, if we know that $x_1$ and $x_2$ are normalized ($\vec{x_1}^T\vec{x_1} = \vec{x_2}^T\vec{x_2} = 1$), then the combined weight is indeed the mean of the individual weights.
A: Consider you would have the following exact relationship between your inputs and your outputs.
$$y = 1 \cdot x_1 + 2 \cdot x_2 + \varepsilon.$$
How, will you be able to fit one single coefficient such that
$$y = \beta (x_1 + x_2) + \varepsilon$$
is still an exact relationship?
You can also proof this by determining the coefficients. If we apply the method of least squares we obtain
$$\hat{\beta}_3 = \dfrac{\overline{yx_1} + \overline{yx_2}}{\overline{x^2_1} + 2 \overline{x_1x_2} + \overline{x^2_2}}$$
Apply the method of least squares to the first regression and you will see that the resulting parameters do not sum up to $\hat{\beta}_3$.
A: Because there is an error term:)
Try:
> y <- x1 + x2 + rnorm(100, 0, 0)           # epsilon = 100 zeroes
> m0 <- lm(y ~ x1 + x2)
> m1 <- lm(y ~ x)
> mean(coef(m0)[2:3]) - coef(m1)[2]
            x 
-2.220446e-16 
> y <- x1 + x2 + rnorm(100, 0, 1)           # epsilon = standard normal
> m0 <- lm(y ~ x1 + x2)
> m1 <- lm(y ~ x)
> mean(coef(m0)[2:3]) - coef(m1)[2]
          x 
-0.01250147 
> y <- x1 + x2 + rnorm(100, 0, 5)           # epsilon = large variance
> m0 <- lm(y ~ x1 + x2)
> m1 <- lm(y ~ x)
> mean(coef(m0)[2:3]) - coef(m1)[2]
         x 
0.08345448 
> y <- x1 + x2 + rnorm(100, 0, 200)         # epsilon = huge variance
> m0 <- lm(y ~ x1 + x2)
> m1 <- lm(y ~ x)
> mean(coef(m0)[2:3]) - coef(m1)[2]
        x 
0.1529181 

