Covariance matrix and Cholesky decomposition

Question

I am simulating a spread option with stochastic volatility using Monte Carlo simulation. I have the positive-definite covariance matrix $$ \rho = \left( \begin{array}{cccc} 1 & \rho_{1,2} & \rho_{1,3} & \rho_{1,4} \\ \rho_{2,1} & 1 & \rho_{2,3} & \rho_{2,4} \\ \rho_{3,1} & \rho_{3,2} & 1 & \rho_{3,4}\\ \rho_{4,1} & \rho_{4,2} & \rho_{4,3} & 1 \end{array} \right) $$

which by definition can be decomposed into the product of two matrixes through Cholesky decomposition in the following way: $$ \rho = L L^T $$ where $T$ indicates the transposed matrix.

In the literature, this factorization renders a system of equations of the form: $$ x_1 = z_1 \\ x_2 = \rho_{1,2} z_1 + \sqrt{1 - \rho_{1,2}^2z_2} \\ x_3 = ... $$

My question is the following. Why do we ignore the matrix $L^T$ when writing the final equations for the random deviates? Why isn't this approach leading to wrong equations, given that we are 'ignoring' part of the system? It seems like we are forgetting half of the problem.

I know this is the correct way to compute random deviates, but I would like to know the reason why this approach works.

Answer 1

I am not sure if I understood your question correctly but I will try to answer it anyway.

If you have a standard normal random vector $z \sim N(\mathbb{0},I_n)$ (where $z,0 \in \mathbb{R}^{n\times1}$ and $I_n \in \mathbb{R}^{n\times n}$ is the identity matrix) and you want to transform it into a multivariate normal $x \sim N(\mu,\Sigma)$ you do it the following way: If $L$ is the Cholesky factorization of $\Sigma$, $\Sigma = LL^T$, then

$$ x = \mu + Lz.$$

Why? Because when you want to see if you really end up with a covariance matrix $\Sigma$ in this construction, it works:

$$\mathbb{E}[(x-\mathbb{E}(x))^T(x-\mathbb{E}(x))] = \mathbb{E}[(Lz)^T(Lz)]=L^T\mathbb{E}[z^Tz]L=L^TI_nL = \Sigma.$$

Or lets think the other way around: If you want to see what matrix you have to multiply $z$ with to create a distribution with covariance matrix $\Sigma$, it turns out to be $L$.

Maybe one more thing: Some textbooks even define the multivariate normal like this: Suppose you have a vector $z$ of independent standard normals and a matrix $L$ such that $LL^T$ is SPD then the distribution of $x=\mu + Lz$ is called multivariate normal with $\mu$ and $LL^T=\Sigma$.