Let $\boldsymbol{X} = (X_1,X_2)^{\rm T}\sim \mathcal{N}_2(\boldsymbol{\mu}, \mathrm{\Sigma})$, where
\begin{eqnarray*}
\boldsymbol{\mu} = (\mu_1, \mu_2)^{\rm T}& = &(\sqrt{\xi_1\xi_2/(\xi_1+\xi_2)}, 0)^{\rm T}\\
\mathrm{\Sigma} & = &\begin{pmatrix} 1 & -\rho\\
-\rho & 1\end{pmatrix}\\
\rho & = & \sqrt{\xi_1\xi_3/(\xi_1+\xi_2)(\xi_2+\xi_3)}.
\end{eqnarray*}
It is given that $\xi_1\leq\xi_2\leq\xi_3$, $\xi_i\geq 0$ and $\sum_{i=1}^3\xi_i = 1$. I have a function
\begin{equation*}
\pi(\boldsymbol{\mu};\boldsymbol{\xi}) = 1-\mathbb{P}(\boldsymbol{X}\leq \boldsymbol{t}) = 1-\mathbb{P}(X_1\leq t \cap X_2\leq t),
\end{equation*}
where $t>1$. Under the constraints $\xi_1\leq\xi_2\leq\xi_3$, $\xi_i\geq 0$ and $\sum_{i=1}^3\xi_i = 1$, numerically we are getting the maximum of $\pi(\boldsymbol{\mu};\boldsymbol{\xi})$ where the non-zero mean $\mu_1$ is maximized. This occurs when $\xi_1 = \xi_2 = \xi_3 = 1/3$.
Can we use some kind of stochastic ordering arguments to prove the result theoretically?
Probability – Problem Related to Bivariate Normal Stochastic Order
pr.probabilityprobability distributions
Related Solutions
$\newcommand{\der}{\mathrm{der}}\newcommand{\tdert}{\mathrm{dert}}\newcommand{\erf}{\operatorname{erf}}\newcommand{\eqs}{\overset{\text{sign}}=}\newcommand{\tder}{\widetilde\der}$The answer is no.
Indeed, let $x:=\xi_1$ and $y:=\xi_2$, so that $\xi_3=1-x-y$, $0\le x\le y\le1-x-y$, whence $x\in[0,1/3]$.
Consider further the case $y=x\in(0,1/3]$, so that $$\mu_1=M(x):=\sqrt{x/2},$$ $$\rho=-R(x),\quad R(x):=\sqrt{\frac{1-2x}{2(1-x)}}.$$ Let also $$p_t(x):=1-P(X_1\le t,X_2\le t).$$ Your desired result would imply that \begin{equation*} p_t(x)\le p_t(1/3) \tag{0} \end{equation*} for all $x\in(0,1/3]$ and all real $t>0$.
Note that \begin{equation*} p_t(x)=1-P(X_1\le t,X_2\le t)=P(M(x),-R(x)), \tag{1} \end{equation*} where \begin{equation*} P(m,r):=1-\int_{-\infty}^t du \int_{-\infty}^t dv\, f_{m,r}(u,v) \end{equation*} and $f_{m,r}$ is the density function of the bivariate normal distribution with means $m,0$, variances $1,1$, and correlation $r$.
The key is Plackett's observation (formula (3)) that \begin{equation*} D_r f_{m,r}(u,v)=D_v D_u f_{m,r}(u,v), \end{equation*} where $D_w$ denote the partial derivative with respect to a variable $w$. It follows that \begin{equation*} \begin{aligned} D_r P(m,r)&:=-\int_{-\infty}^t du \int_{-\infty}^t dv \, D_r f_{m,r}(u,v) \\ &=-\int_{-\infty}^t du \int_{-\infty}^t dv \, D_v D_u f_{m,r}(u,v) \\ &=-f_{m,r}(t,t). \end{aligned} \tag{2} \end{equation*} Next, \begin{equation*} \begin{aligned} P(m,r)&=1-\int_{-\infty}^t du \int_{-\infty}^t dv\, f_{0,r}(u-m,v) \\ &=1-\int_{-\infty}^{t-m} dw \int_{-\infty}^t dv\, f_{0,r}(w,v) \end{aligned} \end{equation*} and hence \begin{equation*} \begin{aligned} D_m P(m,r)=\int_{-\infty}^t dv\, f_{0,r}(t-m,v)=\int_{-\infty}^t dv\, f_{m,r}(t,v); \end{aligned} \tag{3} \end{equation*} the latter integral can be easily expressed in terms of the error function $\erf$ and elementary functions.
By (1) and a chain rule of differentiation, \begin{equation*} D_x p_t(x)=D_m P(M(x),-R(x))M'(x)-D_r P(M(x),-R(x))R'(x). \end{equation*}
In particular, \begin{equation} D_x p_{1/10}(x)\big|_{x=1/3}=\erf\left(\frac{1}{60} \left(3 \sqrt{6}-10\right)\right)+1-\frac{3 }{\sqrt{\pi }}\,e^{(30 \sqrt{6}-77)/1800} =-0.739\ldots<0. \end{equation} So, inequality (0) fails to hold for $t=1/10$ and all $x$ in a left neighborhood of $1/3$. $\quad\Box$
The joint cdf of a multivariate distribution uniquely determines the distribution. So, if the cdf's of multivariate normal distributions are the same, then their mean vectors must be the same (and their covariance matrices must be the same).
To simplify the notation, let $d:=\delta$ and $b_i:=\xi_i$.
Then your three multivariate normal distributions will be the same if their covariance matrices are the same and, for instance, $$d=0,\ b_1=0,\ b_2=\frac9{34},\ b_3=\frac{25}{34}$$ -- with unequal $b_i$'s, contrary to your desired conclusion.
If we impose the condition $d\ne0$, then your three multivariate normal distributions can never be the same.
Here are the calculations in Mathematica:
Best Answer
$\newcommand{\der}{\mathrm{der}}\newcommand{\tdert}{\mathrm{dert}}\newcommand{\erf}{\operatorname{erf}}\newcommand{\eqs}{\overset{\text{sign}}=}\newcommand{\tder}{\widetilde\der}$The statement about the maximum is true, but the proof has hardly anything to do with stochastic ordering arguments. Rather, the problem is reduced to certain problems of real algebraic geometry -- that is, to solving (rather complicated) systems of algebraic inequalities over $\mathbb R$.
Let $x:=\xi_1$ and $y:=\xi_2$, so that $\xi_3=1-x-y$, $0\le x\le y\le1-x-y$, $$\mu_1=M(x,y):=\sqrt{\frac{x y}{x+y}},$$ $$\rho=-R(x,y),\quad R(x,y):=\sqrt{\frac{x(1-x-y)}{(x+y)(1-x)}}.$$ Let also $$p_t(x,y):=1-P(X_1\le t,X_2\le t)$$ and \begin{equation*} G:=\{(x,y)\colon 0\le x\le y\le1-x-y\}=\{(x,y)\colon 0\le x\le1/3,x\le y\le(1-x)/2\}. \tag{-1} \end{equation*} We want to show that \begin{equation*} p_t(x,y)\le p_t(1/3,1/3) \tag{0} \end{equation*} for all $(x,y)\in G$ and all real $t\ge1$.
Note that \begin{equation*} p_t(x,y):=1-P(X_1\le t,X_2\le t)=P(M(x,y),-R(x,y)), \tag{1} \end{equation*} where \begin{equation*} P(m,r):=1-\int_{-\infty}^t du \int_{-\infty}^t dv\, f_{m,r}(u,v) \end{equation*} and $f_{m,r}$ is the density function of the bivariate normal distribution with means $m,0$, variances $1,1$, and correlation $r$.
The key is Plackett's observation (formula (3)) that \begin{equation*} D_r f_{m,r}(u,v)=D_v D_u f_{m,r}(u,v), \end{equation*} where $D_w$ denote the partial derivative with respect to a variable $w$. It follows that \begin{equation*} \begin{aligned} D_r P(m,r)&:=-\int_{-\infty}^t du \int_{-\infty}^t dv \, D_r f_{m,r}(u,v) \\ &=-\int_{-\infty}^t du \int_{-\infty}^t dv \, D_v D_u f_{m,r}(u,v) \\ &=-f_{m,r}(t,t)<0. \end{aligned} \tag{2} \end{equation*} Next, \begin{equation*} \begin{aligned} P(m,r)&=1-\int_{-\infty}^t du \int_{-\infty}^t dv\, f_{0,r}(u-m,v) \\ &=1-\int_{-\infty}^{t-m} dw \int_{-\infty}^t dv\, f_{0,r}(w,v) \end{aligned} \end{equation*} and hence \begin{equation*} \begin{aligned} D_m P(m,r)=\int_{-\infty}^t dv\, f_{0,r}(t-m,v)=\int_{-\infty}^t dv\, f_{m,r}(t,v)>0. \end{aligned} \tag{3} \end{equation*}
Also, for $(x,y)\in G\setminus\{(0,0)\}$, \begin{equation*} (D_x M(x,y),D_y M(x,y))=\frac{(y^2,x^2)}{(x + y)^2}, \end{equation*} \begin{equation*} (D_x R(x,y),D_y R(x,y))=\frac{((1-2 x-y)y,-(1-x) x)}{(1-x)^2 (x + y)^2}, \end{equation*} so that $D_x M(x,y)\ge0$ and $D_x R(x,y)\ge0$.
So, by (1), a chain rule of differentiation, and the inequalities in (2) and (3), \begin{equation*} D_x p_t(x,y)=D_m P(M(x,y),-R(x,y))D_x M(x,y)-D_r P(M(x,y),-R(x,y))D_x R(x,y)\ge0. \end{equation*}
So, the maximum of $p_t(x,y)$ over $(x,y)\in G$ occurs at one of right boundaries of $G$, that is, where either $y=x$ and $y=(1-x)/2$.
Next, for $x\in[0,1/3]$ (cf. (-1)), \begin{equation*} \frac d{dx}\,M(x,(1-x)/2)=\frac{1-2x-x^2}{(1+x)^2}>0, \end{equation*} \begin{equation*} \frac d{dx}\,R(x,(1-x)/2)=\frac1{(1+x)^2}>0. \end{equation*} So, by (1) and the inequalities in (2) and (3), $p_t(x,(1-x)/2)$ is increasing in $x\in[0,1/3]$.
It remains to show that \begin{equation*} \der(t,x)\overset{\text{(?)}}\ge0 \text{ for $x\in[0,1/3]$ and $t\ge1$,} \tag{4} \end{equation*} where \begin{equation*} \begin{aligned} \der(t,x)&:=8 \sqrt{\pi } \sqrt{x}\, e^{\left(\sqrt{2} \sqrt{x}-2 t\right)^2/8}\, D_x p_t(x,x) \\ &=1+\erf\left(v(t,x)/2\right)-\frac{2 }{\sqrt{\pi } (1-x)}\,\sqrt{\frac{x}{1-2 x}}\,e^{u(t,x)/4}, \\ u(t,x)&:=t^2 \left(-4 \sqrt{2} \sqrt{\frac{1-2x}{1-x}} (1-x)+8 x-6\right) \\ &+2 t \left(-2 \sqrt{2} x^{3/2}+2 (1-x) \sqrt{\frac{1-2x}{1-x}} \sqrt{x}+\sqrt{2} \sqrt{x}\right)+x (2 x-1), \\ v(t,x)&:=\sqrt{1-x} \left(\sqrt{\frac{1-2x}{1-x}} \left(\sqrt{2} t-\sqrt{x}\right)+2 t\right). \end{aligned} \end{equation*}
By a chain rule of differentiation, \begin{equation*} D_x p_t(x,x)=D_m P(M(x),-R(x))M'(x)-D_r P(M(x),-R(x))R'(x), \end{equation*} where \begin{equation*} M(x):=M(x,x),\quad R(x):=R(x,x). \end{equation*} Using (2) and (3), we get \begin{equation*} \der(t,x)=A_1(x)e^{A_2(t,x)}+1+\erf(A_3(t,x));\tag{4.5} \end{equation*} here and in what follows, $A_1,A_2,\dots$ are certain algebraic functions and $\erf$ is the error function. So, \begin{equation*} D_t\der(t,x) =A_4(t,x)e^{A_2(t,x)}+A_5(t,x)e^{-A_3(t,x)^2/2}. \tag{5} \end{equation*} thus getting rid of the function $\erf$. Moreover, \begin{equation*} A_5(t,x)>0 \end{equation*} (for $t$ and $x$ as in (4)).
Letting now \begin{equation*} \tdert(t,x):=\frac{D_t\der(t,x)}{A_5(t,x)e^{-A_3(t,x)^2/2}}, \end{equation*} we get \begin{equation*} D_t\der(t,x)\eqs\tdert(t,x)=A_6(t,x)e^{A_7(t,x)}+1, \tag{6} \end{equation*} where $a\eqs b$ means that $a$ and $b$ are of the same sign; this reduces the two exponential expressions $e^{A_2(t,x)}$ and $e^{-A_3(t,x)^2/2}$ in (5) to one such expression. Now we have \begin{multline*} D_t\tdert(t,x)=e^{A_7(t,x)}(1-x) \sqrt{\frac{(1-2x)(1-x)}{x}} \Big(2+\sqrt2\,\sqrt{\frac{1-2x}{1-x}}\,\Big) \\ \times \big(3-4 x+2 \sqrt{2} \sqrt{(1-2 x) (1-x)}\,\big)>0 \end{multline*} (for $x\in(0,1/3]$).
Also,
\begin{equation*} \tdert(0,x) =\frac{(1-2 x) \left(2-2 x+\sqrt{2} \sqrt{(1-2 x) (1-x)}\right)}{\left(2+\sqrt{2} \sqrt{2+\dfrac{1-2 x}{1-x}}\right) (1-x)^2}>0. \end{equation*} So, $\tdert(t,x)>0$ for all $t\ge0$ and all $x\in(0,1/3]$.
So, by (6), it is enough to show that \begin{equation*} \der_1(x):=\der(1,x)\overset{\text{(?)}}\ge0 \tag{7} \end{equation*} for all $x\in(0,1/3]$.
Using (4.5) and some algebra, we get \begin{equation*} \der_1(x)=1-\frac{2/\sqrt\pi}{1-x}\,\sqrt{\frac x{1-2x}}\,e^{p(x)/4}+\erf(q(x)), \end{equation*} where \begin{equation*} p(x):=-6+2 \sqrt{2} \sqrt{x}+7 x-4 \sqrt{2} x^{3/2}+2 x^2+\sqrt{\frac{1-2 x}{1-x}} \left(-4 \sqrt{2}+4 \sqrt{x}+4 \sqrt{2} x-4 x^{3/2}\right), \end{equation*} \begin{equation*} q(x):=\frac{1}{2} \sqrt{1-2 x} \left(\sqrt{2}-\sqrt{x}\right)+\sqrt{1-x}. \end{equation*}
Let \begin{equation*} \der_{11}(x):=\der_1'(x)2\sqrt\pi\, \frac{(1-2 x) (1-x)^2 \sqrt{(1-2x)x} } {p_1(x)e^{p(x)/4}}, \end{equation*} where \begin{equation*} \begin{aligned} p_1(x)&:=(1-2x)(1-x) x p'(x)+2+2x-8 x^2 \\ &=2 + 2 x - 8 x^2 + x (1 - 2 x) ((7 + \sqrt2/\sqrt{x} - 6 \sqrt2 \sqrt{x} + 4 x) (1 - x) +p_2(x)), \end{aligned} \end{equation*} \begin{equation*} p_2(x):=2 (\sqrt2 - \sqrt{x})/ s(x) + (4 \sqrt2 + 2/\sqrt{x} - 6 \sqrt{x}) (1 - x) s(x), \end{equation*} \begin{equation*} s(x):=\sqrt{\frac{1-2 x}{1-x}}. \end{equation*} It is elementary to show that $p_2(x)\ge5$ (for $x\in(0,1/3]$), and hence \begin{equation*} \begin{aligned} p_1(x)&\ge2 + 2 x - 8 x^2 + x (1 - 2 x) ((7 + \sqrt2/\sqrt{x} - 6 \sqrt2 \sqrt{x} + 4 x) (1 - x) +5)>0. \end{aligned} \end{equation*} So, \begin{equation*} \der_{11}(x)\eqs\der_1'(x). \end{equation*} Next, let \begin{equation*} \der_{111}(x):=\der_{11}'(x) \sqrt{(1-2x)x} \, e^{p(x)/4+q(x)^2} \, p_1(x)^2 =\tder_{111}(\sqrt x), \end{equation*} where \begin{equation*} \begin{aligned} \tder_{111}(z)&:= -3 \sqrt{2}+50 z+56 \sqrt{2} z^2-282 z^3-61 \sqrt{2} z^4 \\ &+392 z^5-284 \sqrt{2} z^6+248 z^7+540 \sqrt{2} z^8 \\ &-832 z^9-160 \sqrt{2} z^{10}+384 z^{11}-64 \sqrt{2} z^{12} \\ &-2 \sqrt{\frac{1-2z^2}{1-z^2}} \\ &\times\big(1-6 \sqrt{2} z-32 z^2+48 \sqrt{2} z^3 +33 z^4 \\ & -106 \sqrt{2} z^5 +136 z^6+48\sqrt{2} z^7-230 z^8 \\ &+80 \sqrt{2} z^9+48 z^{10}-64 \sqrt{2} z^{11}+32 z^{12}\big). \end{aligned} \end{equation*} The function $\tder_{111}$ switches its sign only once on the interval $(0,1/\sqrt3]$, from $-$ to $+$, and hence $\der_{11}'$ has the same sign pattern on the interval $(0,1/3]$. So, for some $x_*\in(0,1/3)$, the function $\der_{11}$ is decreasing on $(0,x_*]$ and increasing on $[x_*,1/3]$. Also, $\der_{11}(0+)=-3/2<0$ and $\der_{11}(1/3)=-1.11\ldots<0$. So, $\der_{11}<0$ on $(0,1/3]$. So, $\der_1$ is decreasing on $(0,1/3]$, to $\der_1(0)=1.31\ldots>0$.
Thus, (7) is proved. $\quad\Box$