Simulate random variables with "inner" and "outer" correlation

Question

Let us say we have data grouped in $m$ different classes, each of size $n_j$ for $j = 1,...,m$. We denote as $X_k^{(j)}$ the $k$-th member of group $j$. We want to simulate unit-variance random variables with a structure such that the "inner" correlation in group $j$ is $\rho_j$, and the "outer" correlation between groups $j$ and $i$ is $\rho_{ij}$, with the $\rho$'s as parameters. This is, $$ \begin{aligned} Var(X_k^{(j)}) &= 1 \quad \forall k,j \\ Cov(X_k^{(j)}, X_r^{(j)}) &= \rho_j \quad \forall k \neq r, j = 1,...,m\\ Cov(X_k^{(j)}, X_r^{(i)}) &=\rho_{ij} \quad \forall i\neq j, k\neq r \end{aligned} $$ Basically, I have grouped data and want to create a relationship structure both inside the groups and between the groups. I have tried approaching the problem by using the Gaussian copula and establishing $$ X_k^{(j)} = \sqrt{\beta_j}\, Z + \sqrt{\rho_j - \beta_j} \, Z_j + \sqrt{1-\rho_j} \, e_k^{(j)} $$ assuming both inner and outer correlations are positive, with $Z, Z_j, e_k^{(j)}$ i.i.d. $N(0,1)$. However, one finds that $$ \beta_i \beta_j = \rho_{ij}^2 $$ which I believe is not generally solvable, given there are more equations than variables. I even used non-linear optimization techniques to try and find an approximate solution, without luck (cannot approximate each $\rho_{ij}$ appropriately).

I would like to know if there is any literature regarding this kind of problems or models, or if it is even possible. Even though $\rho_j$ seem to be always positive in my data, I would like the model to be as flexible as possible.

Answer 1

I found an answer, which was easier than expected, using factor models for credit risk. Given $m$ groups each with $n_j$ members, for $j=1,...,m$, we define $$ X_k^{(j)} = \sqrt{\rho_j}Z_j + \sqrt{1-\rho_j}e_k^{(j)} $$ for $k=1,...,n_j$, where $e_k^{(j)}$ are i.i.d. for each and all $k$ and $j$, and independent of $Z_j$ (idiosyncratic noise). Moreover, $\pmb{Z}= (Z_1,...,Z_m)\sim N(0, \Sigma)$ where the entries of $\Sigma$ are given by $$ (\Sigma)_{ij} = \frac{\rho_{ij}}{\sqrt{\rho_i \rho_j}}. $$ Thus, note that $E[X_k^{(j)}] = 0$ and $Var(X_k^{(j)}) = 1$, given the independence of $Z_j$ and $e_k^{(j)}$. Moreover, given the linearity of the normal distribution, $X_k^{(j)} \sim N(0,1)$.

For the correlations, note that (we obviate subindexes for they are not necesarry) $$ \begin{aligned} Corr(X^{(i)}, X^{(j)}) &= Cov(X^{(i)}, X^{(j)}) \quad \text{(Unitary variance)} \\ &= Cov(\sqrt{\rho_i}Z_i + \sqrt{1-\rho_i}e^{(i)},\sqrt{\rho_j}Z_j + \sqrt{1-\rho_j}e^{(j)}) \\ &= \sqrt{\rho_i}\sqrt{\rho_j}Cov(Z_i, Z_j) \quad \text{(Linearity of $Cov$ and ind.)} \\ &= \sqrt{\rho_i}\sqrt{\rho_j}(\Sigma)_{ij} \\ &= \sqrt{\rho_i}\sqrt{\rho_j}\frac{\rho_{ij}}{\sqrt{\rho_i \rho_j}} \\ &= \rho_{ij} \end{aligned} $$ as required. When both the variables are from the same group, we just get $Corr(X^{(j)}, X^{(j)}) = \rho_{jj}$ so we can establish $\rho_{jj} = \rho_j$ w.l.o.g.. In fact, one can prove that $\pmb{X} = (X_1^{(1)},X_2^{(1)},...,X_{n_1}^{(1)},X_1^{(2)},..,X_{n_m}^{(m)})$ the whole vector of $n = \sum_j n_j$ variables is multivariate normal.

Therefore, to simulate such a $\pmb{X}$ one can just construct the whole covariance matrix of size $n$ and simulate a multivariate normal from that (which is cumbersome given that $n$ is large in our application), or you could use the approach above and simulate an $m$ sized vector for $\pmb{Z}$ with covariance matrix $\Omega$, and $n$ i.i.d. $N(0,1)$. We have found that, for our application, the later approach is much faster given that $n$ is in the millions and $m$ around 10.