Sample many correlated random matrices all with the same pairwise correlation coefficient

Question

I am looking to generate $K$ different correlated random matrices, of which the elements all have the same pairwise correlation coefficient. That is, let $A_1, A_2, \ldots, A_K$ be $N \times N$ random matrices with mean $\mu$ and variance $\sigma^2$. Then I want that the the pairwise correlation coefficient of elements abides: $$ \text{cor}[(a_{ij})_k (a_{ij})_\ell] = \rho \quad \forall i, j $$ and for all $k \neq \ell$.

I know that for $K=2$, I can simply generate $A_1 \sim \mathcal{N}(\mu,\sigma^2)$ and then let $A_2 = \rho A_1 + \sqrt{1 - \rho^2} B$, where $B \sim \mathcal{N}(\mu,\sigma^2)$. However, I fail to see how to generalize this to $K$ matrices. My naive implementation simply uses the above formula for two matrices, but just repeat it, i.e. $$ A_{k+1} = \rho A_k + \sqrt{1-\rho^2} B_k $$ where each $B_k$ is i.i.d. with mean $\mu$ and variance $\sigma^2$. But of course, this does not lead to the correct pairwise correlations between all matrices. I also tried Cholesky decomposition with some covariance matrix $\Sigma$, but I fail to see what dimensions $\Sigma$ should have when I want to generate $N\times N$ matrices.

How should I proceed?

Answer 1

Most methods assume your $\mathbf{A}_1, \mathbf{A}_2, \ldots, \mathbf{A}_K$ matrices are instead $n$-length vectors $\mathbf{x}_1, \mathbf{x}_2, \ldots, \mathbf{x}_p$, which form an $n \times p$ data matrix $\mathbf{X}$. Correlation between all possible $(j,k)$ pairs of vectors will result in a $p \times p$ measured correlation matrix $\mathbf{R}$.

Fundamentally, you really can't "say" or state what correlation coefficients should be via a specified (theoretical) multivariate correlation matrix $\boldsymbol{\rho}$ and expect them to appear (i.e., be the same) in a measured (empirical) correlation matrix $\mathbf{R}$. This is because:

1. your data vectors are not infinitely large
2. $\boldsymbol{\rho}$ may violate Cauchy-Schwarz inequality
3. $\boldsymbol{\rho}$ may contain pathologies
4. $\boldsymbol{\rho}$ may be positive semi-definite instead of
   positive definite.

The measured correlation matrix $\mathbf{R}$ does not suffer from any of the above issues, since it's calculated from empirical data vectors $\mathbf{x}_1, \mathbf{x}_2, \ldots, \mathbf{x}_p$.

Since the relationship you are using, i.e., $z=\rho x + \sqrt{1-\rho^2} y$ blows up with more than 2-3 vectors, if you want to continue to use what you have implemented thus far, then you will need to look at methodological procedures published in ref. 1 below (Rebonato et al.) in order to remove pathologies from a specified $\boldsymbol{\rho}$ matrix.

If you want to be able to simulate correlation matrices without pathologies, then use e.g. the Iman & Conover method2. But it still helps to remove pathologies in any simulated correlation matrix. If you want all the correlation matrices to be similar, then just run Iman & Conover's method multiple times using the same specified $\boldsymbol{\rho}$ matrix to generate sets of $\mathbf{x}_j$ which result in sets of measured correlation matrices $\mathbf{R}_1, \mathbf{R}_2, \ldots, \mathbf{R}_K$. Recall, however, none of the resulting $\mathbf{R}_j$ matrices will have the same correlation coefficients.

Summing everything up, when simulating a correlation matrix I either set $\boldsymbol{\rho} \leftarrow \mathbf{R}$ based on measured correlation and use Iman & Conover, or specify $\boldsymbol{\rho}$ and clean it up using methods described in Rebonato et al., then run Iman & Conover.

References:

1. R. Rebonato, P. J¨ackel. The most general methodology to create a valid correlation matrix for risk management and option pricing purposes. J. of Risk. 2:17–27; 2000.

2. R.L. Iman, W.J. Conover. A distribution-free approach to inducing rank correlation among input variables. Commun. Statist. Simulat. Computat. 11(3): 311–334; 1982.