Suppose we have a vector of random variables $X \sim N(\mu_0, \Sigma_0)$. We are interested in the distribution of $Y$ defined as follows:

$Y = (X - \mu_0)^T\Sigma_0^{-1}(X - \mu_0) - (X - \mu_1)^T\Sigma_1^{-1}(X - \mu_1)$

Note that the distribution of $(X - \mu_0)^T\Sigma_0^{-1}(X - \mu_0)$ is Chi-squared distribution. The question is that what is the distribution of $Y$.

Thank you.

  • $\begingroup$ Since the two terms both involve the same $X$, they certainly won't be independent in general, so that's not something we can "note." Precisely what assumptions are you making about the relationship between $\Sigma_0$ and $\Sigma_1$? $\endgroup$
    – whuber
    Commented Sep 17, 2015 at 15:05
  • $\begingroup$ You are right. We can assume that $\Sigma_0$ and $\Sigma_1$ are independent (the same for $\mu_0$ and $\mu_1$). $\endgroup$
    – bassir
    Commented Sep 17, 2015 at 15:11
  • 1
    $\begingroup$ What does "independence" of $\Sigma_0$ and $\Sigma_1$ mean? Are these random matrices in some sense? And if so, how does one guarantee that they are invertible and that at least $\Sigma_0$ is positive semidefinite? $\endgroup$ Commented Sep 17, 2015 at 15:34
  • $\begingroup$ By "independence", I just mean that there is no relationships between $\Sigma_0$ and $\Sigma_1$. $\Sigma_0$ and $\Sigma_1$ are invertible and positive semidefinite. $\endgroup$
    – bassir
    Commented Sep 17, 2015 at 16:07

1 Answer 1


In case you wanted $Y$ to be defined as follows:

$$ X_i \sim N(\mu_i, \Sigma_i) \; \text{for} \;i={0,1} $$ so that $$ Y = (X_0-\mu_0)'\Sigma_0^{-1}(X_0-\mu_0) \; - \; (X_1-\mu_1)'\Sigma_1^{-1}(X_1-\mu_1) $$ where the subtrahend and minuend both are $\sim \chi^2(n)$. I therefore assume that both $X_i$ follow a $n$-variate normale distribution and that $X_0$ and $X_1$ are independent. Then the difference $Y$ follows the variance-gamma distribution. There was already a similar question to this problem in the math community on stackexchange:

Transformation of Difference of Random Variables

In R you can verfiy the result you find following the above link via the following code chunk:


sigma0 <- diag(4)
sigma1 <- diag(4)*2 + 1

X0 <- mvtnorm::rmvnorm(50000, mean=rep(0, 4), sigma0)
X1 <- mvtnorm::rmvnorm(50000, mean=rep(0, 4), sigma1)

X <- cbind(X0, X1)

Y <- apply(X, 1, function(z) { t(z[1:4])%*%solve(sigma0)%*%z[1:4] - t(z[5:8])%*%solve(sigma1)%*%z[5:8] } )

truehist(Y, nbins = 50, ylim=c(0,0.2))
lines(x=seq(-20,20,0.01)-4, y=dvg(x=seq(-20,20,0.01), vgC=4, sigma=4, theta=0, nu=0.5))

enter image description here

  • $\begingroup$ The problem is that in my case, $X_0 = X_1$. $\endgroup$
    – bassir
    Commented Sep 19, 2015 at 13:35
  • $\begingroup$ Let's assume $\Sigma_0 = \Sigma_1$. By this assumption, I think that the problem becomes easier $\endgroup$
    – bassir
    Commented Sep 22, 2015 at 8:10

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.