Vandermonde matrix, $V_n(a_1, \dots, a_n)$ = $\left|\begin{array}{cccc}1 & a_1 & a^2_1 & ... & a^{n-1}_1 \\... & ... & ... & ... \\1 & a_n & a_n^2 & ... & a^{n-1}_n\end{array}\right|$. Let $P(x) = V_n(x,a_2, \dots, a_n)$. $P(x) = (-1)^{n-1}V_{n-1}\prod_{i=2}^n(x-a_i)$.
What I don't see is how to produce the recursion formula $V_n(a_1,\dots, a_n) = \prod_{i=2}^n (a_i - a_1)V_{n-1}(a_2,\dots,a_n)$.
Secondly, how would you use the recursion formula to deduce that $V_n(a_1,\dots, a_n) = \prod_{i,j=1, i < j}^n (a_j-a_i)$?
Do the operations $C_i\leftarrow C_i-a_1C_{i-1}$ for $2\leqslant i\leqslant n$ (this does not change the determinant). The first row is $\pmatrix{1&0&\dots&0}$, and we have \begin{align} V(a_1,\dots,a_n)&=\det\pmatrix{1&0&0&\dots&0\\ 1&a_2-a_1&a_2^2-a_1a_2&\dots& a_2^{n-1}-a_1a_2^{n-2}\\ \vdots& \vdots&&&\vdots\\ 1&a_n-a_1&a_n^2-a_1a_n&\dots&a_n^{n-1}-a_1a_n^{n-1}} \\ &=\det\pmatrix{1&0&0&\dots&0\\ 1&a_2-a_1&a_2(a_2-a_1)&\dots& a_2^{n-2}(a_2-a_1)\\ \vdots& \vdots&&&\vdots\\ 1&a_n-a_1&a_n(a_n-a_1)&\dots&a_n^{n-2}(a_n-a_1)}. \end{align} The formula $V_n(a_1,\dots, a_n) = \prod_{i=2}^n (a_i - a_1)V_{n-1}(a_2,\dots,a_n)$ now follows from an expansion with respect tothe first row and multilinearity of the determinant.