Here is what I found using Mathematica. First, the constant of integration is found for the truncated bivariate distribution, then the first two moments are found followed by the variance. Then the variance is simplified such that we end up with the usual notation.
Here are the results:
$$E(X|Y<y)=\mu_X-\frac{\rho \sigma_X \phi \left(\frac{y-\mu_Y}{\sigma_Y}\right)}{\Phi \left(\frac{y-\mu_Y}{\sigma_Y}\right)}$$
$$V(X|Y<y)=\sigma_X^2 \left(\frac{\rho ^2 \phi \left(\frac{y-\mu_Y}{\sigma_Y}\right) \left((\mu_Y-y) \Phi \left(\frac{y-\mu_Y}{\sigma_Y}\right)-\sigma_Y \phi \left(\frac{y-\mu_Y}{\sigma_Y}\right)\right)}{\sigma_Y \Phi \left(\frac{y-\mu_Y}{\sigma_Y}\right)^2}+1\right)$$
If $\mu_X=\mu_Y=1$ and $\sigma_X=\sigma_Y=1$, then
$$E(X|Y<y)=-\frac{\rho \phi (y)}{\Phi (y)}$$
$$V(X|Y<y)=1-\frac{\rho ^2 \phi (y) (y \Phi (y)+\phi (y))}{\Phi (y)^2}$$
Here is the code:
(* Define bivariate normal distribution *)
d = BinormalDistribution[{μX, μY}, {σX, σY}, ρ];
(* Get constant of integration for truncated distribution *)
c = 1/Integrate[PDF[d, {x, y}], {y, -∞, y0}, {x, -∞, ∞},
Assumptions -> y0 ∈ Reals && σX > 0 && σY > 0 && μX ∈ Reals && μY ∈ Reals && -1 < ρ < 1];
(* Mean *)
ex = c Integrate[x PDF[d, {x, y}], {y, -∞, y0}, {x, -∞, ∞},
Assumptions -> y0 ∈ Reals && σX > 0 && σY > 0 && μX ∈ Reals && μY ∈ Reals && -1 < ρ < 1];
(* Expectation of X^2 *)
ex2 = c Integrate[x^2 PDF[d, {x, y}], {y, -∞, y0}, {x, -∞, ∞},
Assumptions -> y0 ∈ Reals && σX > 0 && σY > 0 && μX ∈ Reals && μY ∈ Reals && -1 < ρ < 1];
(* V(X|Y<y *)
var = ex2 - ex^2
(* Now attempt to simplify and write in terms of usual notation *)
var = var // FullSimplify
var = var /. Erfc[z_] -> 1 - Erf[z] //. Erf[Abs[z_]/(Sqrt[2] σY)] Sign[z_]^3 ->
Erf[z/(Sqrt[2] σY)]
var = var /. E^(-((y0 - μY)^2/(2 σY^2))) -> Sqrt[2 π] ϕ[(y0 - μY)/σY] /.
Erf[z_] -> -1 + 2 Φ[Sqrt[2] z] /. y0 -> y // FullSimplify
(* E(X|Y<y) *)
expectation = ex /. E^(-((y0 - μY)^2/(2 σY^2))) -> Sqrt[2 π] ϕ[(y0 - μY)/σY] /.
Erf[z_] -> -1 + 2 Φ[Sqrt[2] z] /. y0 -> y // FullSimplify