It's more than odd. The two approaches you give should give the same result but only the second one below does give the right answer. (You should inform Wolfram, Inc.)
FullSimplify[Probability[x y >= 0, {x, y} \[Distributed] MultinormalDistribution[{{1, \[Rho]}, {\[Rho], 1}}]],
Assumptions -> -1 <= \[Rho] <= 1]
(* 1/2 *)
2 Probability[x >= 0 && y >= 0, {x, y} \[Distributed] MultinormalDistribution[{{1, \[Rho]}, {\[Rho], 1}}]] // FullSimplify
(* 1/2 + ArcSin[\[Rho]]/\[Pi] *)
An alternative way that does work with x*y is the following:
distxy = TransformedDistribution[x y, {x, y} \[Distributed] MultinormalDistribution[{{1, \[Rho]}, {\[Rho], 1}}]]]
1 - CDF[distxy, 0] // FullSimplify
(* 1/2 + ArcSin[\[Rho]]/\[Pi] *)
