The Danielson-Lanczos lemma is the basis for fast Fourier transform algorithms.
Now, I do understand this step
$\displaystyle X_{k} = \sum_{n=0}^{N-1} x_{n}\omega^{kn}_{N} = \sum_{n=0}^{(N/2)-1} x_{2n}\omega^{k(2n)}_{N} + \sum_{n=0}^{(N/2)-1} x_{2n+1}\omega^{k(2n+1)}_{N}$
but I don't see how the rest follows.
So, keeping that $\omega_N = e^{-\frac{2\pi i}{N}}$, $\hat{x}_n = x_{2n}$, $\tilde{x}_n = x_{2n+1} \; (n = 0, 1, \dots, \frac{N}{2}-1)$, let's regroup terms further:
$$ \begin{align} & X_{k} = \sum_{n=0}^{N-1} x_{n}\omega^{kn}_{N} = \sum_{n=0}^{(N/2)-1} x_{2n}\omega^{k(2n)}_{N} + \sum_{n=0}^{(N/2)-1} x_{2n+1}\omega^{k(2n+1)}_{N} = \\ & \sum_{n=0}^{(N/2)-1} \hat{x}_{n}\omega^{2kn}_{N} + \sum_{n=0}^{(N/2)-1} \tilde{x}_{n}\omega^{2kn}_{N}\cdot\omega^k_N = \sum_{n=0}^{(N/2)-1} \hat{x}_{n}\omega^{kn}_{N/2} + \omega^k_N \cdot \sum_{n=0}^{(N/2)-1} \tilde{x}_{n}\omega^{kn}_{N/2} = \\ & \hat{X}_k + \omega^k_N \cdot \tilde{X}_k , \end{align}$$ where $\hat{X}$ and $\tilde{X}$ denote Fourier transforms of $\hat{x}$ and $\tilde{x}$ (they're both $N/2$-length vectors). My $\hat{x}$ and $\tilde{x}$ stand for $F^e$ and $F^o$ in description from your link.