joint distribution, probability, calculating probabilities under false independent assumption, when the random variables are actually dependent

Suppose I have random variables $X,Y,Z$ and I would like to compute the probability that random variable $X$ is smaller than $Y$ and $Z$:

$$\pi_X \overset{def}{=} Pr(X < Y, X < Z) = \int Pr(x < Y,x < Z)dPr(X\le x)$$ Define $\pi_Y$ and $\pi_Z$ similarly and without loss of generality, let $\pi_X \le \pi_Y \le \pi_Z$. I know that $X,Y$ and $Z$ are positively associated. Suppose that I compute $\pi_X$ "falsely" assuming that $Z \perp Y$ using the formula below:

$$\pi^{ind}_X = \int Pr(x < Y)Pr(x < Z)dPr(X\le x).$$

I would like to know what condition I need for the the dependence structure among $X,Y$ and $Z$ (e.g., positive correlation?) to claim $\pi_X \ge \pi_X^{ind}$? I.e., the $\pi_X^{ind}$ under the (false) independence assumption always underestimates the chance that $X$ is smaller than both $Y$ and $Z$.

Based on the following small simulation study in R, I am guessing that $\pi_X \ge \pi_X^{ind}$ is true when $X$, $Y$ and $Z$ are "positively" associated in some way. The simulation model assume that $X$, $Y$ and $Z$ are from the multivariate normal as:

$$\begin{bmatrix} X\\ Y\\ Z \end{bmatrix} \sim N_3\left( \begin{bmatrix} -2\\ 0\\ 2 \end{bmatrix} \sim \begin{bmatrix} \sigma^2&\rho\sigma &\rho\sigma\\ \rho\sigma& 1 &\rho\\ \rho\sigma & \rho & 1\\ \end{bmatrix} \right)$$ where I consider the values $\sigma=0.5,1,1.5,2$. Then I generate large enough samples from this multivariate distribution and computed $\pi_X$ based on the Monte Carlo samples. I also computed $\pi_X^{ind}$ assuming the independence across $X,Y$ and $Z$. I see that regardless of the values of $\sigma$, the $\pi_X$ is always greater than $\pi_X^{ind}$ when $rho > 0$.

library(mvtnorm)
get.pi <- function(xyz)
{
## obtain the MC estimate of pi^ind and pi
xs <- xyz[,1]
ys <- xyz[,2]
zs <- xyz[,3]
num.indp <- num.dp <- 0
M <- length(xs)
for (ix in 1 : M)
{
x.temp <- xs[ix]
num.dp <- num.dp + mean(x.temp < ys & x.temp < zs)
num.indp <- num.indp + mean(x.temp < ys)*mean(x.temp < zs)
}
return(list(indp=num.indp/M,dp=num.dp/M))
}

N <- 10000
sds1 <- seq(0.5,2,0.5)
sd2 <- 1
sd3 <- 1
rhos <- seq(0,1,0.05)
mus <- c(-2,0,2)
cols <- c("red","green")
par(mfrow=c(2,length(sds1)/2),mar=c(2,2,2,2))
for (isd in 1 : length(sds1))
{
sd1 <-sds1[isd]
res.indP <- res.noIndP <- rep(NA, length(rhos))
for (irho in 1 : length(rhos))
{
rho <- rhos[irho]
rho12 <- rho
rho13 <- rho
rho23 <- rho

S <- matrix(c(sd1^2,sd1*sd2*rho12,sd1*sd3*rho13,
sd1*sd2*rho12,sd2^2,sd2*sd3*rho13,
sd1*sd3*rho13,sd2*sd3*rho13,sd3^2),3,3)

xyz <- rmvnorm(N,mean=mus,sigma=S)
re <- get.pi(xyz)
res.indP[irho] <- re$indp res.noIndP[irho] <- re$dp
}
rang <- c(res.indP,res.noIndP)
ylim <- c(min(rang),max(rang))
plot(rhos,res.indP,type="l",col=cols[1],ylim=ylim,
main=paste("sd1=",sd1,sep=""))
points(rhos,res.noIndP,type="l",col=cols[2],ylim=ylim)
}

legend("topright",c("Independent","No assumption"),col=cols,lwd=1)
• I believe the answer will vary depending on the precise meaning of "positively associated." What exactly do you mean by this? – whuber Nov 14 '14 at 3:18
• @whuber I am actually looking which condition of "positive association" is required to satisfy $\pi_X \ge \pi_X^{ind}$... – FairyOnIce Nov 14 '14 at 3:20
• You might want to emphasize that by editing your question. – whuber Nov 14 '14 at 3:22
• @whuber Thank you for the advice. I just edited my question. – FairyOnIce Nov 14 '14 at 3:25
• This post in math.SE, math.stackexchange.com/questions/452484/…, may be interesting as a more general treatment of such issues. – Alecos Papadopoulos Nov 14 '14 at 10:23