Given any $n-1$-dimensional linear subspace $U \subseteq \mathbb{R}^n$ orthogonal to a unit vector $p$, how to generate Gaussian-distributed vectors on $U$ with covariance $\Sigma \in \mathbb{R}^{(n-1)} \times \mathbb{R}^{(n-1)}$, centered at the origin?

Something like this in the 3D case:

enter image description here

  • 1
    $\begingroup$ are you asking for the math or code? $\endgroup$ – Haitao Du May 18 '17 at 16:59
  • $\begingroup$ How exactly is the covariance given to you? Is it an $n\times n$ semi-definite singular matrix? If so, what do you perceive as being the difference between this situation and one where $\Sigma$ is an $n\times n$ definite nonsingular matrix? (Many methods for the latter, usual case work just fine in the former case.) $\endgroup$ – whuber May 18 '17 at 17:01
  • $\begingroup$ I am asking for a general procedure/pseudo code. The covariance is a non-singular $(n-1)\times(n-1)$ matrix. $\endgroup$ – Haochi Kiang May 18 '17 at 17:04
  • $\begingroup$ My point is that virtually any procedure that works for nonsingular matrices will work, or can be made to work, with singular ones. If your favorite procedure fails, add a tiny positive multiple of the identity matrix to $\Sigma$: that will not appreciably change the distribution but it will be positive-definite and nonsingular. $\endgroup$ – whuber May 18 '17 at 17:09
  • $\begingroup$ But since $\Sigma$ is $n-1$ dimensional, say, if I do mvrnorm(1, rep(0,n-1), Sigma) in R, I will end up with a $n-1$-vector, not a $n$-vector. I tried adding a 0 to the nth dimension to get a $n$-vector and then project it onto $U$, but this doesn't work for all $p$. $\endgroup$ – Haochi Kiang May 18 '17 at 17:17
  1. Diagonalize $\Sigma$, i.e. find its eigenvalues $\sigma_i^2$ (which will be non-negative) and eigenvectors $\vec{v}_i$.
  2. Generate $(n-1)$ normal deviates $z_i \sim N(0, \sigma_i)$. Form $\vec{u} = \sum_{i} z_i \vec{v}_i$. This gives you a multi-normally distributed vector $u$ in the subspace.
  3. If necessary, express the subspace vector $u$ in terms of the full space coordinates. This is just vector rotation and addition.

The details depend on whether $\Sigma$ is expressed as an $(n-1) \times (n-1)$ matrix in the subspace coordinate system or whether $\Sigma$ is expressed as an $n \times n$ matrix in the full coordinate system. In the former case, all eigenvalues will be positive and you will need step (3). In the latter case, one eigenvalues will be zero, the others will span the subspace, and you will not need step (3).

  • $\begingroup$ "Just vector rotation and addition" in (3), which is the essential part of this answer, belies an important lacuna in the question: that information is not present in an $n-1\times n-1$ covariance matrix. It has to be supplied in addition to it. $\endgroup$ – whuber May 18 '17 at 18:07

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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