Hi, I am trying to compare two methods for the calculation of a trivariate normal cdf: 1) "direct" - using the original covariance matrix 2) "indirect" - transforming the integration interval using the Mahalanobis transformation and then integrating as a standard random variable. The following code demonstrates the two methods. Unfortunately, I get different results. Am I missing something? Thanks in advance. Shlomi Lifshits Tel-Aviv University #################### library(mvtnorm) #create a symmetric positive definite matrix for the covariance matrix my_mat<-matrix(data=c(1,2,3,0,0.1,0,0,0,0.1), nrow=3, ncol=3, byrow=TRUE) covmat<-my_mat%*%t(my_mat) #calculate the Mahalanobis transformation spectral_dec<-eigen(covmat) new_eigenval<-c(sqrt(spectral_dec$values[1]),sqrt(spectral_dec$values[2]),sqrt(spectral_dec$values[3])) covmat_trans<-solve(spectral_dec$vectors%*%diag(new_eigenval)%*%t(spectral_dec$vectors)) # integrate with the original limit orig_limit<-c(-0.5,0.5,0.5) a<-pmvnorm(lower=-Inf, upper=orig_limit, mean=rep(0, 3), sigma=covmat) #transform the limit and integrate new_limit_raw<-covmat_trans%*%orig_limit new_limit<-c(new_limit_raw[1,1],new_limit_raw[2,1],new_limit_raw[3,1]) new_sigma<-diag(3) b<-pmvnorm(lower=-Inf, upper=new_limit, mean=rep(0, 3), sigma=new_sigma) # a,b should be the same but they are not!!!?