Here is the definition of joint normality in my textbook.
Def: Two random variables $X$ and $Y$ are said to be jointly normal if they have the joint density $$f_{X,Y}(x,y)\\=\frac{1}{2\pi\sigma_1\sigma_2\sqrt{1-\rho^2}}\exp\left\{-\frac{1}{2(1-\rho^2)}\left[ \frac{(x-\mu_1)^2}{\sigma_1^2} - \frac{2\rho(x-\mu_1)(y-\mu_2)}{\sigma_1 \sigma_2} +\frac{(y-\mu_2)^2}{\sigma_2^2} \right] \right\},$$
where $\sigma_1 >0, \sigma_2>0,|\rho|<1,$ and $\mu_1,\mu_2$ are real numbers.
My textbook states without proof that this definition is equivalent to the statement that linear combinations of jointly normal random variables are still jointly normal.
I tried to google a proof for this equivalence but I didn't manage to find one. Can anyone help me with proving the $(\Longrightarrow)$ direction? Thanks.
Textbook page:
Stochastic calculus for finance II Continuous time models, Steven E. Shreve.