Show that $E(Y\mid X=x)$ is a linear function in $x$ 

Let $Y$ and $X$ be bivariate normal distributed with expectationvector $\mu=(\mu_Y,\mu_X)^T$ and covariance matrix $\Sigma=\begin{pmatrix}\sigma_Y^2 & p_{XY}\\p_{XY} & \sigma_X^2\end{pmatrix}$. Show that the conditional expectation $E(Y\mid X=x)$ is a linear function in $x$.


Hello!
To my knowledge it is
$$
E(Y\mid X=x)=\int_{\mathbb{R}}y\cdot f_{Y\mid X}(y\mid x)\, dy.
$$
So first I tried to determine $f_{Y\mid X}(y\mid x)$ by
$$
f_{Y\mid X}(y\mid x)=\frac{f_{Y,X}(y,x)}{f_X(x)}.
$$
To my calculation this is
$$
f_{Y\mid X}(y\mid x)=\frac{\sigma_X}{\sqrt{2\pi}\sqrt{\sigma_Y^2\sigma_X^2-p_{XY}^2}}\exp\left(-\frac{1}{2(\sigma_Y^2\sigma_X^2-p_{XY}^2)}\cdot(\sigma_X^2(y-\mu_Y)^2-2p_{XY}(x-\mu_X)(y-\mu_Y)+\sigma_Y^2(x-\mu_X)^2)+\frac{1}{2}\frac{(x-\mu_X)^2}{\sigma_X^2}\right)
$$
Is that right?
In case it is: How can I know determine
$$
\int y\cdot f_{Y\mid X}(y\mid x)\, dy,
$$
i.e. how can I determine
$$
\frac{\sigma_X}{\sqrt{2\pi}\sqrt{\sigma_Y^2\sigma_X^2-p_{XY}^2}}\int_{\mathbb{R}}y\cdot\exp\left(-\frac{1}{2(\sigma_Y^2\sigma_X^2-p_{XY}^2)}\cdot(\sigma_X^2(y-\mu_Y)^2-2p_{XY}(x-\mu_X)(y-\mu_Y)+\sigma_Y^2(x-\mu_X)^2)+\frac{1}{2}\frac{(x-\mu_X)^2}{\sigma_X^2}\right)\, dy?
$$
Edit:
If I set 
$$
c:=-\frac{\sigma_X^2}{2(\sigma_Y^2\sigma_X^2-p_{XY}^2)}, d:=\frac{2p_{XY}(x-\mu_X)}{2(\sigma_Y^2\sigma_X^2-p_{XY}^2)}, q:=\frac{\sigma_Y^2(x-\mu_X)^2}{2(\sigma_Y^2\sigma_X^2-p_{XY}^2)}
$$
and
$$
w:=\frac{(x-\mu_X)^2}{2\sigma_X^2}
$$
then I have to calculate the following:
$$
\frac{\sigma_X\exp(-q+w)}{\sqrt{2\pi}\sqrt{\sigma_Y^2\sigma_X^2-p_{XY}^2}}\int y\cdot\exp(c(y-\mu_Y)^2+d(y-\mu_Y))\, dy
$$
 A: You might want to work directly on random variables: there exists some parameters $(a,b)$ and a standard normal random variable $Z$ independent of $X$ such that $Y-\mu_Y=a(X-\mu_X)+bZ$ (can you show this?), thus $E(Z\mid X)=E(Z)=0$ hence $E(Y\mid X)=\mu_Y+a(X-\mu_X)$. To compute $a$, note that $\mathrm{Cov}(Y,X)=a\cdot\mathrm{var}(X)+b\cdot\mathrm{Cov}(Z,X)$ and $\mathrm{Cov}(Z,X)=0$, hence $\rho_{XY}=a\cdot\sigma^2_X$.
A: If you insist on working with densities, then note that the joint density of $(X,Y)$ is given by
$$
f_{X,Y}(x,y)=\frac{1}{2\pi\sigma_X\sigma_Y\sqrt{1-\rho^2}}\exp\Big(-\frac{1}{2(1-\rho^2)}r(x,y)\Big),
$$
with
$$
r(x,y)=\frac{(x-\mu_X)^2}{\sigma_X^2}-2\rho\frac{(x-\mu_X)}{\sigma_X}\frac{(y-\mu_Y)}{\sigma_Y}+\frac{(y-\mu_Y)^2}{\sigma_Y^2}.
$$
Here $\rho$ is the correlation between $X$ and $Y$, i.e. $p_{XY}=\rho \sigma_X\sigma_Y$. Calculations (do these) now show that $r(x,y)$ can be written as
$$
r(x,y)=-\frac{1}{2\sigma_X^2}(x-\mu_X)^2-\frac{1}{2\sigma_Y^2(1-\rho^2)}\Big(y-\mu_Y-\frac{\rho\sigma_Y}{\sigma_X}(x-\mu_X)\Big)^2
$$
and hence $f_{X,Y}$ factors into the product of
$$
\frac{1}{\sqrt{2\pi}\sigma_X}\exp\Big(-\frac{1}{\sigma_X^2}(x-\mu_X)^2\Big)\tag{1}
$$
and
$$
\frac{1}{\sqrt{2\pi}\sigma_Y(1-\rho^2)}\exp\left(-\frac{1}{2\sigma_Y^2(1-\rho^2)}\Big(y-\mu_Y-\frac{\rho\sigma_Y}{\sigma_X}(x-\mu_X)\Big)^2\right)\tag{2}
$$
showing that the conditional density of $Y$ given $X=x$ is given by $(2)$, i.e. $Y\mid X=x$ is normally distributed with mean $\mu_Y+\rho\frac{\sigma_Y}{\sigma_X}(x-\mu_X)$ and variance $\sigma_Y^2(1-\rho^2)$.
A: If you want to calculate directly, try massaging your expression for 
$f_{Y\mid X}(y\mid X=x)$ into the form
$$f_{Y\mid X}(y\mid X=x) = \frac{\exp\left(-\frac{1}{2}\left(\frac{y-a}{b}\right)^2\right)}{b\sqrt{2\pi}}$$
where $a$ and $b$ are functions of $x$ and the other parameters but
do not depend on $y$. Then see if you can determine the expected
value without doing the integration that you propose to carry out.
