# Find a design matrix giving a prescribed covariance matrix

Given a symmetrical positive-definite matrix $L$, how do I find a matrix $X$ such that $L = X^T X$?

• This isn't really on-topic, but the Cholesky decomposition is typically the most convenient. May 29 '15 at 8:33
• It is not possible to recover design matrix $X$ from the product $X^TX$, without any additional information, since for any orthogonal matrix $O$, the product $(OX)^T(OX)=X^TX$ does not change. May 29 '15 at 9:59
• @Dougal What is the reason you think this question is not on-topic here?
– whuber
May 29 '15 at 13:37
• @whuber I somehow missed the title of the question in the review page, and thought it was purely a linear algebra question. With the edited title it's clearer that it's appropriate here. May 29 '15 at 17:03

Essentially, you can reverse the principal components process.

Let's start by making some uncorrelated, centered and normalized data, you can create the columns in your matrix by drawing from any distribution as long as the draws are independent. Here I draw from a centered and standardized uniform distribution $X \sim U(-.5, .5) / \sqrt{12}$, but only because it is convenient and illustrates the process:

> X <- matrix(runif(120000, -.5, .5)/ (1/sqrt(12)),
nrow=40000)
> t(X) %*% X / 400000
[,1]          [,2]          [,3]
[1,]  0.0996386570 -0.0002036308 -0.0006817355
[2,] -0.0002036308  0.1000222332  0.0015659048
[3,] -0.0006817355  0.0015659048  0.1000916918


That is:

$$\frac{1}{N} X^t X = I$$

Now let's write down a symmetric, positive definite matrix that we want to be our final correlation matrix:

S <- matrix(c(  1, -.5,   0,
-.5,   1, -.5,
0, -.5,   1), nrow=3, byrow=TRUE
)


Use the sqrtm function to compute the matrix square root:

> rtS <- sqrtm(S)
> t(rtS) %*% rtS
[,1] [,2]          [,3]
[1,]  1.000000e+00 -0.5  5.551115e-17
[2,] -5.000000e-01  1.0 -5.000000e-01
[3,]  5.551115e-17 -0.5  1.000000e+00


Then form the matrix product $X\sqrt{S}$. This has the property that:

$$(X\sqrt{S})^t X\sqrt{S} = \sqrt{S}^t X^t X \sqrt{S} = N \sqrt{S}^t \sqrt{S} = NS$$

so

$$\frac{1}{N} (X\sqrt{S})^t (X\sqrt{S}) = S$$

so $X \sqrt{S}$ is the data matrix you want

> t(X %*% sqrtm(S)) %*% (X %*% sqrtm(S)) / 40000
[,1]       [,2]         [,3]
[1,]  0.998555833 -0.4987640 -0.009608847
[2,] -0.498764047  0.9921826 -0.483641573
[3,] -0.009608847 -0.4836416  0.993164671