Como gerar uma matriz de correlação aleatória grande com várias correlações fortes presentes?

25

Gostaria de gerar uma matriz de correlação aleatório de n x n de tamanho de tal modo que existem algumas correlações moderadamente fortes presente:Cn×n

  • matriz quadrada de reais simétrico de tamanho, por exemplo, com n = 100 ;n×nn=100
  • positivo-definido, ou seja, com todos os autovalores reais e positivos;
  • classificação completa;
  • todos os elementos diagonais são iguais a ;1 1
  • elementos fora da diagonal devem ser razoavelmente distribuídos uniformemente em . A distribuição exata não importa, mas eu gostaria de ter uma quantidade moderadamente grande (por exemplo, 10 % ) de valores moderadamente grandes (por exemplo, com valor absoluto de 0,5 ou superior). Basicamente eu quero ter certeza de que C é não quase diagonal com todos os elementos fora da diagonal 0 .(-1 1,1 1)10%0,5C0 0

Existe uma maneira simples de fazer isso?

O objetivo é usar essas matrizes aleatórias para comparar alguns algoritmos que trabalham com matrizes de correlação (ou covariância).


Métodos que não funcionam

Aqui estão algumas maneiras de gerar matrizes de correlação aleatória que eu conheço, mas que não funcionam para mim aqui:

  1. Gere aleatório de tamanho s × n , centralize, padronize e forme a matriz de correlação C = 1Xs×n. Ses>n, isso geralmente resultará em todas as correlações fora da diagonal em torno de0. Ses«n, algumas correlações será forte, masCnão será posto completo.C=1s1XXs>n0snC

  2. Gere a matriz definida positiva aleatória de uma das seguintes maneiras:B

    • Gere o quadrado aleatório e faça B positivo definido simétrico B = A A positive .AB=AA

    • Gere o quadrado aleatório , torne simétrico E = A + A e torne-o positivo definitivo, realizando a decomposição do Eigen E = U S U e definindo todos os autovalores negativos para zero: B = UAE=A+AE=USU . NB: isso resultará em uma matriz deficiente na classificação.B=Umax{S,0}U

    • Gere ortogonal aleatório (por exemplo, gerando quadrado aleatório A e fazendo sua decomposição QR, ou via processo de Gram-Schmidt) e diagonal D aleatória com todos os elementos positivos; forma B = Q D Q .QADB=QDQ

    Obteve matriz pode ser facilmente normalizada ter todos aqueles na diagonal: C = D - 1 / 2 B D - 1 / 2 , onde D = d i um gBC=D1/2BD1/2 é a matriz diagonal com o mesmo diagonal como B . Todas as três maneiras listadas acima para gerar B resultam em C com elementos fora da diagonal perto de 0 .D=diagBBBC0


Atualização: Tópicos mais antigos

Depois de postar minha pergunta, encontrei duas quase duplicatas no passado:

Infelizmente, nenhum desses tópicos continha uma resposta satisfatória (até agora :)

ameba diz Restabelecer Monica
fonte
11
Você pode criar matriz ortogonal aleatória por processos QR ou Gram-Schmidt. Isso será "autovetores do PCA". Adicione escala às suas colunas (transforme em "cargas"). Obtenha a matriz de covariância desses carregamentos. Algo assim ...
ttnphns
11
Uhm, bem .. Imagine que queremos criar uma nXkmatriz de carregamento W, não totalmente aleatória, mas a que queremos (ela WW'+diag(noise)definirá a matriz cov que procuramos. A única tarefa é corrigir o W normalizado pela coluna (ou seja, k "eigenvectors") para se tornar ortogonal Qualquer método para de-correlato variáveis correlacionadas (aqui variáveis são os autovetores) provavelmente irá fazer (Esta uma idéia crua)...
ttnphns
11
Ah, @whuber, agora entendo o que você quer dizer. Sim, você está certo: se todos os elementos fora da diagonal forem idênticos e iguais a , então a matriz é de fato a posição completa e com definição positiva ... Isso não é, obviamente, o que eu tinha em mente: eu gostaria da distribuição de elementos fora da diagonal em cada matriz devem ser razoavelmente "espalhados", não a distribuição entre matrizes ...ρ
ameba diz Reinstate Monica
3
Você pode querer olhar para a distribuição LKJ
shadowtalker
2
@ttnphns: Acho que finalmente entendi que você estava certa o tempo todo: o que você sugeriu é a maneira mais simples de atingir a meta. Eu adicionei uma atualização à minha resposta implementando essencialmente o que você escreveu acima.
Ameba diz Reinstate Monica

Respostas:

14

Outras respostas surgiram com bons truques para resolver meu problema de várias maneiras. No entanto, encontrei uma abordagem baseada em princípios que acredito ter uma grande vantagem de ser conceitualmente muito clara e fácil de ajustar.

Nesta discussão: Como gerar eficientemente matrizes de correlação aleatória positiva-semidefinida? - Descrevi e forneci o código para dois algoritmos eficientes de geração de matrizes de correlação aleatória. Ambos vêm de um artigo de Lewandowski, Kurowicka e Joe (2009), ao qual o @ssdecontrol se refere nos comentários acima (muito obrigado!).

Por favor, veja minha resposta lá para muitas figuras, explicações e código do matlab. O chamado método "vine" permite gerar matrizes de correlação aleatórias com qualquer distribuição de correlações parciais e pode ser usado para gerar matrizes de correlação com grandes valores fora da diagonal. Aqui está a figura de exemplo desse segmento:

Método Vine

A única coisa que muda entre as subparcelas é um parâmetro que controla quanto a distribuição das correlações parciais se concentra em torno de .±1 1

Copio meu código para gerar essas matrizes aqui também, para mostrar que não é mais longo que os outros métodos sugeridos aqui. Por favor, veja minha resposta vinculada para algumas explicações. Os valores de betaparampara a figura acima foram (e a dimensionalidadefoi 100 ).50.,20,10,5,2,1 1d100

function S = vineBeta(d, betaparam)
    P = zeros(d);           %// storing partial correlations
    S = eye(d);

    for k = 1:d-1
        for i = k+1:d
            P(k,i) = betarnd(betaparam,betaparam); %// sampling from beta
            P(k,i) = (P(k,i)-0.5)*2;     %// linearly shifting to [-1, 1]
            p = P(k,i);
            for l = (k-1):-1:1 %// converting partial correlation to raw correlation
                p = p * sqrt((1-P(l,i)^2)*(1-P(l,k)^2)) + P(l,i)*P(l,k);
            end
            S(k,i) = p;
            S(i,k) = p;
        end
    end

    %// permuting the variables to make the distribution permutation-invariant
    permutation = randperm(d);
    S = S(permutation, permutation);
end

Atualização: autovalores

@psarka pergunta sobre os autovalores dessas matrizes. Na figura abaixo, planto os espectros de autovalor das mesmas seis matrizes de correlação acima. Observe que eles diminuem gradualmente; por outro lado, o método sugerido por @psarka geralmente resulta em uma matriz de correlação com um grande autovalor, mas o restante é bastante uniforme.

autovalores das matrizes acima


Atualizar. Método realmente simples: vários fatores

Semelhante ao que @ttnphns escreveu nos comentários acima e @GottfriedHelms em sua resposta, uma maneira muito simples de alcançar meu objetivo é gerar aleatoriamente várias cargas de fator ( )k<n (matriz aleatória detamanho k × n ), formar o matriz de covariância W W (que obviamente não será a classificação completa) e adicione a ela umamatriz diagonal aleatóriaWk×nWW com os elementos positivos para fazer B = W W + DDB=WW+Dclassificação completa. A matriz de covariância resultante pode ser normalizada para se tornar uma matriz de correlação (conforme descrito na minha pergunta). Isso é muito simples e faz o truque. Aqui estão alguns exemplos de matrizes de correlação para :k=100,50.,20,10,5,1 1

matrizes de correlação aleatória a partir de fatores aleatórios

A única desvantagem é que a matriz resultante terá valores próprios grandes e, em seguida, uma queda repentina, em oposição a uma boa decadência mostrada acima com o método vine. Aqui estão os espectros correspondentes:k

eigenspectra dessas matrizes

Aqui está o código:

d = 100;    %// number of dimensions
k = 5;      %// number of factors

W = randn(d,k);
S = W*W' + diag(rand(1,d));
S = diag(1./sqrt(diag(S))) * S * diag(1./sqrt(diag(S)));
ameba diz Restabelecer Monica
fonte
+1. No entanto, aqui está apenas um lembrete para sua última seção sobre "método fatorial". A abordagem estritamente correta chama que as colunas de Wsão ortogonais (isto é, cossenos entre elas são 0). Simplesmente gerar aleatoriamente, Wclaro, não fornece isso. Se eles não são ortogonais - ou seja, são factores oblíqua (chamada então Wcomo W_) - factor de teorema não é WW'mas W_CW_'com Csendo "correlações" (co-senos) entre factores. Agora, C=Q'Qcom Qsendo a matriz de rotação não ortogonal de rotação W_=inv(Q)'W(e assim W=W_Q'). Gere um pouco Q- uma matriz com a coluna ss = 1 e matriz ss = tamanho da matriz.
ttnphns
... erro de digitação: não W_=inv(Q)'W, é claro W_= W inv(Q)'.
ttnphns
@ttnphns: O que você está dizendo está correto, mas não acho que isso importe com o objetivo de gerar matrizes de correlação aleatória. Se eu gerar aleatoriamente, sim, suas colunas não serão exatamente ortogonais, mas W WWWW+DW
11
Traduzindo isso para R:W = replicate(k, rnorm(d)); S = W%*%t(W) + diag(rnorm(d),nrow=d); S = diag(1/sqrt(diag(S)))%*%S%*%diag(1/sqrt(diag(S)))
Scott Worland
11
@ Mihai, bom ponto e suas sugestões são provavelmente as mais simples. Você também pode fazerS <- matrix(nearPD(S, corr = TRUE, keepDiag = TRUE)$mat@x,ncol(S),ncol(S))
Scott Worland
7

uma

import numpy as np
from random import choice
import matplotlib.pyplot as plt

n = 100
a = 2

A = np.matrix([np.random.randn(n) + np.random.randn(1)*a for i in range(n)])
A = A*np.transpose(A)
D_half = np.diag(np.diag(A)**(-0.5))
C = D_half*A*D_half

vals = list(np.array(C.ravel())[0])
plt.hist(vals, range=(-1,1))
plt.show()
plt.imshow(C, interpolation=None)
plt.show()

A distribuição um tanto uniforme Os resultados do imshow

psarka
fonte
crsk[-uma,uma]X
Sim, você está completamente certo! (Oh garoto, isso foi realmente bobo: D). Mudei a parte aleatória para randn (1) * a e agora está muito melhor.
psarka
k
uman
Uma desvantagem desse método é que a matriz de correlação resultante tem um valor próprio grande, mas os demais são quase uniformes. Portanto, esse procedimento não produz uma matriz de correlação "geral" ... Não que eu o tenha especificado na minha pergunta. Mas o @ssdecontrol mencionado nos comentários acima indica que aparentemente existem maneiras de amostrar todas as matrizes de correlação; isso parece interessante, mas muito mais complicado.
Ameba diz Reinstate Monica
6

Hmm, depois de fazer um exemplo na minha linguagem MatMate, vejo que já existe uma resposta em python, que pode ser preferível porque o python é amplamente usado. Mas como você ainda tinha perguntas, mostro a minha abordagem usando a linguagem Matmate-matrix, talvez seja mais auto-comentada.

Método 1
(usando o MatMate):

v=12         // 12 variables
f=3          // subset-correlation based on 3 common factors
vg = v / f   // variables per subsets

 // generate hidden factor-matrix
             // randomu(rows,cols ,lowbound, ubound) gives uniform random matrix 
             //    without explicite bounds the default is: randomu(rows,cols,0,100)
L = {   randomu(vg,f)     || randomu(vg,f)/100  || randomu(vg,f)/100 , _
        randomu(vg,f)/100 || randomu(vg,f)      || randomu(vg,f)/100 , _
        randomu(vg,f)/100 || randomu(vg,f)/100  || randomu(vg,f)     }

 // make sure there is itemspecific variance
 // by appending a diagonal-matrix with random positive entries
L = L || mkdiag(randomu(v,1,10,20)) 
  // make covariance and correlation matrix
cov = L *'   // L multiplied  with its transpose
cor = covtocorr(cov)
                   set ccdezweite=3 ccfeldweite=8
                   list cor
cor = 
   1.000,   0.321,   0.919,   0.489,   0.025,   0.019,   0.019,   0.030,   0.025,   0.017,   0.014,   0.014
   0.321,   1.000,   0.540,   0.923,   0.016,   0.015,   0.012,   0.030,   0.033,   0.016,   0.012,   0.015
   0.919,   0.540,   1.000,   0.679,   0.018,   0.014,   0.012,   0.029,   0.028,   0.014,   0.012,   0.012
   0.489,   0.923,   0.679,   1.000,   0.025,   0.022,   0.020,   0.040,   0.031,   0.014,   0.011,   0.014
   0.025,   0.016,   0.018,   0.025,   1.000,   0.815,   0.909,   0.758,   0.038,   0.012,   0.018,   0.014
   0.019,   0.015,   0.014,   0.022,   0.815,   1.000,   0.943,   0.884,   0.035,   0.012,   0.014,   0.012
   0.019,   0.012,   0.012,   0.020,   0.909,   0.943,   1.000,   0.831,   0.036,   0.013,   0.015,   0.010
   0.030,   0.030,   0.029,   0.040,   0.758,   0.884,   0.831,   1.000,   0.041,   0.017,   0.022,   0.020
   0.025,   0.033,   0.028,   0.031,   0.038,   0.035,   0.036,   0.041,   1.000,   0.831,   0.868,   0.780
   0.017,   0.016,   0.014,   0.014,   0.012,   0.012,   0.013,   0.017,   0.831,   1.000,   0.876,   0.848
   0.014,   0.012,   0.012,   0.011,   0.018,   0.014,   0.015,   0.022,   0.868,   0.876,   1.000,   0.904
   0.014,   0.015,   0.012,   0.014,   0.014,   0.012,   0.010,   0.020,   0.780,   0.848,   0.904,   1.000

O problema aqui pode ser que definimos blocos de submatrizes que possuem altas correlações com pouca correlação entre e isso não é programaticamente, mas pelas constantes expressões de concatenação. Talvez essa abordagem possa ser modelada de forma mais elegante em python.


Método 2 (a)
Depois disso, existe uma abordagem completamente diferente, na qual preenchemos a possível covariância restante por quantidades aleatórias de 100% em uma matriz de cargas fatoriais. Isso é feito no Pari / GP:

{L = matrix(8,8);  \\ generate an empty factor-loadings-matrix
for(r=1,8, 
   rv=1.0;    \\ remaining variance for variable is 1.0
   for(c=1,8,
        pv=if(c<8,random(100)/100.0,1.0); \\ define randomly part of remaining variance
        cv= pv * rv;  \\ compute current partial variance
        rv = rv - cv;     \\ compute the now remaining variance
        sg = (-1)^(random(100) % 2) ;  \\ also introduce randomly +- signs
        L[r,c] = sg*sqrt(cv) ;  \\ compute factor loading as signed sqrt of cv
       )
     );}

cor = L * L~

e a matriz de correlação produzida é

     1.000  -0.7111  -0.08648   -0.7806   0.8394  -0.7674   0.6812    0.2765
   -0.7111    1.000   0.06073    0.7485  -0.7550   0.8052  -0.8273   0.05863
  -0.08648  0.06073     1.000    0.5146  -0.1614   0.1459  -0.4760  -0.01800
   -0.7806   0.7485    0.5146     1.000  -0.8274   0.7644  -0.9373  -0.06388
    0.8394  -0.7550   -0.1614   -0.8274    1.000  -0.5823   0.8065   -0.1929
   -0.7674   0.8052    0.1459    0.7644  -0.5823    1.000  -0.7261   -0.4822
    0.6812  -0.8273   -0.4760   -0.9373   0.8065  -0.7261    1.000   -0.1526
    0.2765  0.05863  -0.01800  -0.06388  -0.1929  -0.4822  -0.1526     1.000

Possivelmente, isso gera uma matriz de correlação com componentes principais dominantes devido à regra de geração cumulativa para a matriz de cargas fatoriais. Também pode ser melhor garantir uma definição positiva, tornando a última parte da variação um fator único. Deixei no programa para manter o foco no princípio geral.

Uma matriz de correlação 100x100 tinha as seguintes frequências de correlações (arredondadas para 1 dec)

    e    f            e: entry(rounded) f: frequency
  -----------------------------------------------------
  -1.000, 108.000
  -0.900, 460.000
  -0.800, 582.000
  -0.700, 604.000
  -0.600, 548.000
  -0.500, 540.000
  -0.400, 506.000
  -0.300, 482.000
  -0.200, 488.000
  -0.100, 464.000
   0.000, 434.000
   0.100, 486.000
   0.200, 454.000
   0.300, 468.000
   0.400, 462.000
   0.500, 618.000
   0.600, 556.000
   0.700, 586.000
   0.800, 536.000
   0.900, 420.000
   1.000, 198.000

[atualizar]. Hmm, a matriz 100x100 está mal condicionada; Pari / GP não pode determinar os valores próprios corretamente com as polroots (charpoly ()) - função mesmo com precisão de 200 dígitos. Fiz uma rotação de Jacobi para formar pca na matriz de loadings L e encontro principalmente autovalores extremamente pequenos, imprimi-os em logaritmos na base 10 (que fornecem aproximadamente a posição do ponto decimal). Leia da esquerda para a direita e depois linha por linha:

log_10(eigenvalues):
   1.684,   1.444,   1.029,   0.818,   0.455,   0.241,   0.117,  -0.423,  -0.664,  -1.040
  -1.647,  -1.799,  -1.959,  -2.298,  -2.729,  -3.059,  -3.497,  -3.833,  -4.014,  -4.467
  -4.992,  -5.396,  -5.511,  -6.366,  -6.615,  -6.834,  -7.535,  -8.138,  -8.263,  -8.766
  -9.082,  -9.482,  -9.940, -10.167, -10.566, -11.110, -11.434, -11.788, -12.079, -12.722
 -13.122, -13.322, -13.444, -13.933, -14.390, -14.614, -15.070, -15.334, -15.904, -16.278
 -16.396, -16.708, -17.022, -17.746, -18.090, -18.358, -18.617, -18.903, -19.186, -19.476
 -19.661, -19.764, -20.342, -20.648, -20.805, -20.922, -21.394, -21.740, -21.991, -22.291
 -22.792, -23.184, -23.680, -24.100, -24.222, -24.631, -24.979, -25.161, -25.282, -26.211
 -27.181, -27.626, -27.861, -28.054, -28.266, -28.369, -29.074, -29.329, -29.539, -29.689
 -30.216, -30.784, -31.269, -31.760, -32.218, -32.446, -32.785, -33.003, -33.448, -34.318

[atualização 2]
Método 2 (b)
Uma melhoria pode ser aumentar a variação específica de itens para algum nível não marginal e reduzir para um número razoavelmente menor de fatores comuns (por exemplo, número inteiro inteiro do número de item):

{  dimr = 100;
   dimc = sqrtint(dimr);        \\ 10 common factors
   L = matrix(dimr,dimr+dimc);  \\ loadings matrix 
                                \\     with dimr itemspecific and 
                                \\          dimc common factors
   for(r=1,dim, 
         vr=1.0;                \\ complete variance per item 
         vu=0.05+random(100)/1000.0;   \\ random variance +0.05
                                       \\ for itemspecific variance
         L[r,r]=sqrt(vu);              \\ itemspecific factor loading  
         vr=vr-vu;
         for(c=1,dimc,
                cv=if(c<dimc,random(100)/100,1.0)*vr;
                vr=vr-cv;
                L[r,dimr+c]=(-1)^(random(100) % 2)*sqrt(cv)
             )
        );}

   cov=L*L~
   cp=charpoly(cov)   \\ does not work even with 200 digits precision
   pr=polroots(cp)    \\ spurious negative and complex eigenvalues...

A estrutura do resultado

em termos de distribuição de correlações:imagem

permanece semelhante (também a descompossibilidade desagradável do PariGP), mas os valores próprios, quando encontrados pela rotação jacobi da matriz de loadings, agora têm uma estrutura melhor.

log_10(eigenvalues):
   1.677,   1.326,   1.063,   0.754,   0.415,   0.116,  -0.262,  -0.516,  -0.587,  -0.783
  -0.835,  -0.844,  -0.851,  -0.854,  -0.858,  -0.862,  -0.862,  -0.868,  -0.872,  -0.873
  -0.878,  -0.882,  -0.884,  -0.890,  -0.895,  -0.896,  -0.896,  -0.898,  -0.902,  -0.904
  -0.904,  -0.909,  -0.911,  -0.914,  -0.920,  -0.923,  -0.925,  -0.927,  -0.931,  -0.935
  -0.939,  -0.939,  -0.943,  -0.948,  -0.951,  -0.955,  -0.956,  -0.960,  -0.967,  -0.969
  -0.973,  -0.981,  -0.986,  -0.989,  -0.997,  -1.003,  -1.005,  -1.011,  -1.014,  -1.019
  -1.022,  -1.024,  -1.031,  -1.038,  -1.040,  -1.048,  -1.051,  -1.061,  -1.064,  -1.068
  -1.070,  -1.074,  -1.092,  -1.092,  -1.108,  -1.113,  -1.120,  -1.134,  -1.139,  -1.147
  -1.150,  -1.155,  -1.158,  -1.166,  -1.171,  -1.175,  -1.184,  -1.184,  -1.192,  -1.196
  -1.200,  -1.220,  -1.237,  -1.245,  -1.252,  -1.262,  -1.269,  -1.282,  -1.287,  -1.290
Elmos de Gottfried
fonte
Muito obrigado! Muito interessante, mas levarei algum tempo para digerir ...
ameba diz Reinstate Monica
Ainda tenho que analisar cuidadosamente sua resposta, mas, enquanto isso, leio um artigo sobre amostragem de matrizes de correlação aleatória, e um dos métodos a partir daí pode ser usado para fazer exatamente o que eu preciso. Eu postei uma resposta aqui, você pode estar interessado em dar uma olhada! Ele está relacionado a uma resposta muito mais detalhada que escrevi em outro tópico.
Ameba diz Reinstate Monica
@amoeba: feliz que você encontrou algo muito bom para você! É uma pergunta interessante, voltarei mais tarde, talvez melhore / adapte os procedimentos do MatMate (e as torne sub-rotinas) de acordo com o artigo em que você trabalhou.
Gottfried Helms
2

Pergunta interessante (como sempre!). Que tal encontrar um conjunto de matrizes de exemplo que exibam as propriedades que você deseja e, em seguida, fazer combinações convexas, pois seUMA e B são definitivos positivos, então também λUMA+(1 1-λ)B. Como bônus, não será necessário redimensionar as diagonais pela convexidade da operação. Ajustando oλpara estar mais concentrado em 0 e 1 em relação à distribuição uniforme, você pode concentrar as amostras nas bordas do politopo ou no interior. (Você pode usar uma distribuição beta / Dirichlet para controlar a concentração versus uniformidade).

Por exemplo, você pode deixar UMA ser simétrico por componentes e Bseja toeplitz. Claro, você sempre pode adicionar outra classeC, e pegue λUMAUMA+λBB+λCC de tal modo que λ=1 1 e λ0 0, e assim por diante.

Andrew M
fonte
Obrigado pela sugestão, Andrew, mas é claro que seria melhor ter um método imparcial que não precise começar com algumas definições predefinidas. UMA e B... Nos comentários ao meu @ssdecontrol pergunta original se refere a um artigo descrevendo algoritmos para matrizes de correlação de amostra uniformemente (em certo sentido preciso), ou inclinado para matriz de identidade, mas não consigo encontrar uma maneira ainda para provar-los inclinado longe de Identidade ... Eu também encontrei alguns tópicos antigos aqui fazendo quase a mesma pergunta, talvez você esteja interessado, veja minha atualização.
Ameba diz Reinstate Monica
Ah, mas a partir de um algoritmo desse tipo e de uma diversidade adequada nos "vértices" (ou seja, matrizes) que definem seu polítopo de matrizes de correlação definida positiva, você pode usar a amostragem por rejeição para obter qualquer distribuição de valores próprios, uniformidade de entradas, etc, que você deseja. No entanto, não está claro para mim qual seria uma boa base. Soa como uma pergunta para alguém que estudou álgebra abstrata mais recentemente do que eu.
Andrew M
Olá novamente, li um artigo sobre amostragem de matrizes de correlação aleatória e um dos métodos a partir daí pode ser usado para fazer exatamente o que eu preciso. Eu postei uma resposta aqui, você pode estar interessado em dar uma olhada! Ele está vinculado a uma resposta muito mais detalhada que escrevi em outro tópico.
Ameba diz Reinstate Monica
2

R possui um pacote (clusterGeneration) que implementa o método em:

Exemplo:

> (cormat10 = clusterGeneration::rcorrmatrix(10, alphad = 1/100000000000000))
        [,1]   [,2]    [,3]     [,4]     [,5]   [,6]   [,7]    [,8]     [,9]   [,10]
 [1,]  1.000  0.344 -0.1406 -0.65786 -0.19411  0.246  0.688 -0.6146  0.36971 -0.1052
 [2,]  0.344  1.000 -0.4256 -0.35512  0.15973  0.192  0.340 -0.4907 -0.30539 -0.6104
 [3,] -0.141 -0.426  1.0000  0.01775 -0.61507 -0.485 -0.273  0.3492 -0.30284  0.1647
 [4,] -0.658 -0.355  0.0178  1.00000  0.00528 -0.335 -0.124  0.5256 -0.00583 -0.0737
 [5,] -0.194  0.160 -0.6151  0.00528  1.00000  0.273 -0.350 -0.0785  0.08285  0.0985
 [6,]  0.246  0.192 -0.4847 -0.33531  0.27342  1.000  0.278 -0.2220 -0.11010  0.0720
 [7,]  0.688  0.340 -0.2734 -0.12363 -0.34972  0.278  1.000 -0.6409  0.40314 -0.2800
 [8,] -0.615 -0.491  0.3492  0.52557 -0.07852 -0.222 -0.641  1.0000 -0.50796  0.1461
 [9,]  0.370 -0.305 -0.3028 -0.00583  0.08285 -0.110  0.403 -0.5080  1.00000  0.3219
[10,] -0.105 -0.610  0.1647 -0.07373  0.09847  0.072 -0.280  0.1461  0.32185  1.0000
> cormat10[lower.tri(cormat10)] %>% psych::describe()
   vars  n  mean   sd median trimmed mad   min  max range skew kurtosis   se
X1    1 45 -0.07 0.35  -0.08   -0.07 0.4 -0.66 0.69  1.35 0.03       -1 0.05

Infelizmente, não parece possível simular correlações que seguem uma distribuição uniforme-ish com isso. Parece fazer correlações mais fortes quando alphadé definido com valores muito pequenos, mas mesmo assim 1/100000000000000, o intervalo de correlações seria de apenas 1,40.

No entanto, espero que isso possa ser útil para alguém.

Deleet
fonte