Stack Exchange

sim.R

library(MASS)

# Generate simulated data
set.seed(1)
mu = c(1,1.5)
Sigma = matrix(c(1.5,-.5,-.5,.7),2,2)
X = mvrnorm(10000000,mu,Sigma)

# Calculate pairwise minimums
mins = pmin(X[,1],X[,2])

# Set minumum
m = .5

# Create subsets for conditional statements
dex = mins < m + .0001 & mins > m - .0001 # min(x,y) = m 
dexx = X[,1] < m + .0001 & X[,1] > m - .0001 # x = m
dexy = X[,2] < m + .0001 & X[,2] > m - .0001 # y = m

# Calculate z values
muxy = mu[1]+Sigma[1,2]*(m-mu[2])/Sigma[2,2]
muyx = mu[2]+Sigma[1,2]*(m-mu[1])/Sigma[1,1]

s2xy = Sigma[1,1] - Sigma[1,2]^2/Sigma[2,2]
s2yx = Sigma[2,2] - Sigma[1,2]^2/Sigma[1,1]

zxy = (m-muxy)/sqrt(s2xy)
zyx = (m-muyx)/sqrt(s2yx)

# Check (1)
mean(X[dex,1]<X[dex,2]) 
mean(X[dexx,2]>m)/(mean(X[dexx,2]>m)+mean(X[dexy,1]>m))
# Check (7)
(1-pnorm(zxy))/(2-pnorm(zxy)-pnorm(zyx))