`ey = resid(lm(y ~ x2 + x3))` calculates the variance of `y` *not explained* by the regressors $x_2$ and $x_3$. On the other hand, `ex = resid(lm(x ~ x2 + x3))` gives you the variance of $x$ *not explained* by $x_2$ and $x_3$.

Therefore, regressing `ey` over `ex` will calculate the 

$\text{variance of }y 
\text{ not explained by }x_2\text{ and} x_3$ 


**explained by** the 

$\text{variance of }x_1 \text{ not explained by }x_2\text{ and} x_3$

So you have eliminated the contribution of $x_2$ and $x_3$, and you are really just calculating the variance of $y$ explained by $x_1$, after having regressed both variables (dependent an independent) over $x_2$ and $x_2$ to effectively eliminate the effect of these two additional explanatory variables. 

Now you just have to keep in mind that the [regression through the origin](http://rinterested.github.io/statistics/OLS_no_linear_algebra.html) with one single regressor is $\hat \beta=\frac{\displaystyle\sum_{i=1}^{n}Y_iX_i}{\displaystyle\sum_{i=1}^{n}X_i^2}$, and that the $-1$ in `coef(lm(ey ~ ex - 1))` is the call for the OLS without intercept.

The reason why the intercept is eliminated is because the intercept is a regressor in its own right. Even if we didn't specify it in the call for `ey` and `ex` the regression over $1$ was there in the model matrix.

The package `{swirl}` in R contains step-by-step practice in this technique of picking one regressor, replacing all other variables by the residuals of their regressions against that one.