TAOCP 3.4.1 Exercise 13
Let $X=(X_1,\ldots,X_n)^T,$ where the $X_i$ are independent normal deviates with mean $0$ and variance $1$.
Section 3.4.1: Numerical Distributions
Exercise 13. [HM25] If a set of $n$ independent normal deviates $X_1, X_2, \ldots, X_n$, with mean 0 and variance 1, show how to find constants $b_i$ and $a_{ij}$, $1 \le i \le n$, so that if
$$Y_1 = b_1 + a_{11}X_1, \quad Y_2 = b_2 + a_{21}X_1 + a_{22}X_2, \quad \ldots, \quad Y_n = b_n + a_{n1}X_1 + \cdots + a_{nn}X_n,$$
then $Y_1, Y_2, \ldots, Y_n$ are dependent normally distributed variables; $Y_i$ has mean $\mu_i$, and the $Y_i$'s have a given covariance matrix $(c_{ij})$. (The covariance, $c_{ij}$, of $Y_i$ and $Y_j$ is defined to be the average value of $(Y_i - \mu_i)(Y_j - \mu_j)$. In particular, $c_{ii}$ is the variance of $Y_i$, the square of its standard deviation. Not all matrices $(c_{ij})$ can be covariance matrices, and your construction is, of course, only supposed to work whenever a solution to the given conditions is possible.)
Verified: yes
Solve time: 1m41s
Solution
Let
$X=(X_1,\ldots,X_n)^T,$
where the $X_i$ are independent normal deviates with mean $0$ and variance $1$. Then
$E(X)=0,\qquad \operatorname{Cov}(X)=I,$
where $I$ is the identity matrix.
We seek constants $b_i$ and $a_{ij}$ such that
$Y=(Y_1,\ldots,Y_n)^T=\mu+AX,$
where
$\mu=(\mu_1,\ldots,\mu_n)^T,$
and
$$ A= \begin{pmatrix} a_{11} & 0 & \cdots & 0\ a_{21} & a_{22} & \cdots & 0\ \vdots & \vdots & \ddots & \vdots\ a_{n1} & a_{n2} & \cdots & a_{nn} \end{pmatrix}. $$
The required representation in the statement of the exercise is exactly the coordinate form of this matrix equation.
Since each $Y_i$ is a linear combination of jointly normal variables, the vector $Y$ is jointly normal. Therefore it suffices to arrange the prescribed means and covariances.
The mean vector is
$$ E(Y)=E(\mu+AX) =\mu+A,E(X) =\mu. $$
Hence we take
$b_i=\mu_i\qquad(1\le i\le n).$
The covariance matrix of $Y$ is
$$ \operatorname{Cov}(Y) = E\bigl((Y-\mu)(Y-\mu)^T\bigr) = E\bigl((AX)(AX)^T\bigr) = A,E(XX^T),A^T. $$
Since $\operatorname{Cov}(X)=I$,
$$ E(XX^T)=I, $$
and therefore
$$ \operatorname{Cov}(Y)=AA^T. $$
Thus the problem reduces to finding a lower triangular matrix $A$ satisfying
$$ AA^T=(c_{ij}). \tag{1} $$
A covariance matrix must be symmetric and nonnegative definite. Conversely, every symmetric nonnegative definite matrix admits a factorization of the form (1). This is the Cholesky factorization.
The coefficients $a_{ij}$ may be computed recursively. Let
$C=(c_{ij}).$
For $i=1,\ldots,n$:
$$ a_{ii} = \sqrt{ c_{ii} - \sum_{k=1}^{i-1}a_{ik}^2 }, \tag{2} $$
and for $j<i$,
$$ a_{ij} = \frac{ c_{ij} - \sum_{k=1}^{j-1}a_{ik}a_{jk} }{ a_{jj} }. \tag{3} $$
These formulas are obtained by equating the entries of $AA^T$ with those of $C$.
Indeed, for $i=j$,
$$ c_{ii} = \sum_{k=1}^{i}a_{ik}^2, $$
which yields (2); and for $i>j$,
$$ c_{ij} = \sum_{k=1}^{j}a_{ik}a_{jk}, $$
which yields (3).
The construction succeeds exactly when $C$ is a valid covariance matrix, because in that case the quantities under the square roots in (2) are nonnegative and a factorization $C=AA^T$ exists.
With $b_i=\mu_i$ and $a_{ij}$ determined by (2) and (3), we obtain
$$ E(Y_i)=\mu_i, $$
and
$$ \operatorname{Cov}(Y_i,Y_j)=c_{ij} \qquad(1\le i,j\le n). $$
Therefore $Y_1,\ldots,Y_n$ are jointly normal random variables having the prescribed means and covariance matrix.
This completes the proof.
∎