Wie wird eine große zufällige Korrelationsmatrix mit vollem Rang und einigen starken Korrelationen erstellt?

25

Ich möchte eine zufällige Korrelationsmatrix C einer Größe von erzeugen, n×nso dass einige mäßig starke Korrelationen vorliegen:

  • quadratische reelle symmetrische Matrix von n×n Größe, mit zB n=100 ;
  • positiv-definit, dh mit allen Eigenwerten real und positiv;
  • voller Rang;
  • alle diagonalen Elemente sind gleich 1 ;
  • Nicht diagonale Elemente sollten einigermaßen gleichmäßig auf (1,1) . Die genaue Verteilung spielt keine Rolle, aber ich hätte gerne eine mäßig große Menge (z. B. 10% ) von mäßig großen Werten (z. B. mit einem absoluten Wert von 0.5 oder höher). Grundsätzlich möchte ich sicherstellen, dass C ist nicht fast diagonal mit allen off-diagonalen Elemente 0 .

Gibt es eine einfache Möglichkeit, dies zu tun?

Der Zweck besteht darin, solche Zufallsmatrizen zu verwenden, um einige Algorithmen zu vergleichen, die mit Korrelations- (oder Kovarianz-) Matrizen arbeiten.


Methoden, die nicht funktionieren

Hier sind einige Möglichkeiten, um zufällige Korrelationsmatrizen zu generieren, von denen ich weiß, die hier jedoch nicht funktionieren:

  1. Erzeuge zufälliges X von s×n Größe, zentriere, standardisiere und bilde die Korrelationsmatrix C=1s1XX. Wenns>n, führt dies im Allgemeinen dazu, dass alle nicht diagonalen Korrelationen um0. Wennsn, sind einige Korrelationen stark, aberChat nicht den vollen Rang.

  2. Generiere zufällige positive definite MatrixB auf eine der folgenden Arten:

    • Erzeuge ein zufälliges Quadrat und mache einen symmetrischen positiven Punkt B = A A .AB=AA

    • Generiere das zufällige Quadrat , mache symmetrisch E = A + A + und mache es positiv, indem du die Eigenzerlegung E = U S U ⊤ durchführst und alle negativen Eigenwerte auf Null setzt: B = UAE=A+AE=USU . Hinweis: Dies führt zu einer Matrix mit Rangdefizit.B=Umax{S,0}U

    • Generiere zufälliges orthogonales (z. B. durch Generieren des zufälligen Quadrats A und Ausführen seiner QR-Zerlegung oder über den Gram-Schmidt-Prozess) und zufällige Diagonale D mit allen positiven Elementen; Form B = Q D Q .QADB=QDQ

    Erhaltenen Matrix leicht normalisiert werden kann alle diejenigen haben auf der Diagonale: C = D - 1 / 2 B D - 1 / 2 , wobei D = d i a gBC=D1/2BD1/2 ist die Diagonalmatrix mit den gleichen Diagonalenwie B . Alle drei oben aufgeführten Wege zur Erzeugung von B führen dazu, dass C nicht diagonale Elemente nahe 0 hat .D=diagBBBC0


Update: Ältere Threads

Nach dem Posten meiner Frage habe ich in der Vergangenheit zwei fast doppelte gefunden:

Leider enthielt keiner dieser Threads eine zufriedenstellende Antwort (bis jetzt :)

Amöbe sagt Reinstate Monica
quelle
1
Sie können eine zufällige orthogonale Matrix durch QR- oder Gram-Schmidt-Prozesse erstellen. Das werden "Eigenvektoren von PCA" sein. Fügen Sie den Spalten eine Skala hinzu (verwandeln Sie sich in "Ladungen"). Holen Sie sich die Kovarianzmatrix aus diesen Ladungen. So etwas in der Art ...
ttnphns 18.11.14
1
Stellen Sie sich vor, wir wollen eine Ladematrix nXkW erzeugen , die nicht völlig zufällig ist, sondern diejenige, die wir wollen (sie wird WW'+diag(noise)die gesuchte cov-Matrix definieren. Die einzige Aufgabe besteht darin, das spaltennormierte W (dh das k) zu korrigieren "Eigenvektoren"), um orthogonal zu werden. Jede Methode zur Dekorrelation korrelierter Variablen (hier sind Variablen die Eigenvektoren) wird dies wahrscheinlich tun. (Dies ist eine
ungefähre
1
Ah, @whuber, jetzt verstehe ich was du meinst. Ja, Sie haben Recht: Wenn alle nicht diagonalen Elemente identisch und gleich , dann ist die Matrix in der Tat voller Rang und positiv-definit ... Dies ist natürlich nicht das, was ich gedacht habe: Ich möchte die Verteilung von Off-Diagonal-Elemente in jeder Matrix werden vernünftigerweise "verteilt", nicht die Verteilung über Matrizen ...ρ
Amöbe sagt Reinstate Monica
3
Vielleicht möchten Sie einen Blick in die LKJ-Distribution
werfen
2
@ttnphns: Ich glaube, ich habe endlich verstanden, dass Sie die ganze Zeit Recht hatten: Was Sie vorgeschlagen haben, ist der einfachste Weg, um zum Ziel zu gelangen. Ich habe meiner Antwort ein Update hinzugefügt, das im Wesentlichen das implementiert, was Sie oben geschrieben haben.
Amöbe sagt Reinstate Monica

Antworten:

14

Andere Antworten lieferten nette Tricks, um mein Problem auf verschiedene Arten zu lösen. Ich habe jedoch einen prinzipiellen Ansatz gefunden, der meiner Meinung nach den großen Vorteil hat, dass er konzeptionell sehr klar und einfach anzupassen ist.

In diesem Thread: Wie können zufällige positiv-semidefinite Korrelationsmatrizen effizient generiert werden? - Ich habe den Code für zwei effiziente Algorithmen zur Erzeugung zufälliger Korrelationsmatrizen beschrieben und bereitgestellt. Beide stammen aus einem Artikel von Lewandowski, Kurowicka und Joe (2009), auf den @ssdecontrol in den obigen Kommentaren Bezug genommen hat (vielen Dank!).

In meiner Antwort finden Sie viele Abbildungen, Erklärungen und Matlab-Code. Das sogenannte "vine" -Verfahren ermöglicht die Erzeugung zufälliger Korrelationsmatrizen mit beliebiger Verteilung von Teilkorrelationen und kann zur Erzeugung von Korrelationsmatrizen mit großen Werten außerhalb der Diagonale verwendet werden. Hier ist das Beispiel aus diesem Thread:

Rebe-Methode

Das einzige, was sich zwischen Unterplots ändert, ist ein Parameter, der steuert, wie stark sich die Verteilung der Teilkorrelationen um .±1

Ich kopiere meinen Code, um diese Matrizen auch hier zu generieren, um zu zeigen, dass er nicht länger ist als die anderen hier vorgeschlagenen Methoden. In meiner verknüpften Antwort finden Sie einige Erklärungen. Die Werte betaparamfür die obige Abbildung betrugen (und die Dimensionalität betrug 100 ).50,20,10,5,2,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

Update: Eigenwerte

@psarka fragt nach den Eigenwerten dieser Matrizen. In der Abbildung unten zeichne ich die Eigenwertspektren derselben sechs Korrelationsmatrizen wie oben. Beachten Sie, dass sie allmählich abnehmen; Im Gegensatz dazu führt die von @psarka vorgeschlagene Methode im Allgemeinen zu einer Korrelationsmatrix mit einem großen Eigenwert, der Rest ist jedoch ziemlich einheitlich.

Eigenwerte der obigen Matrizen


Aktualisieren. Wirklich einfache Methode: mehrere Faktoren

Ähnlich wie @ttnphns in den obigen Kommentaren und @GottfriedHelms in seiner Antwort schrieb, besteht ein sehr einfacher Weg, um mein Ziel zu erreichen, darin, zufällig mehrere ( ) Faktorladungen W (Zufallsmatrix von k × n Größe) zu erzeugen , aus der Kovarianzmatrix W W (die natürlich nicht den vollen Rang hat) und füge dazu eine zufällige Diagonalmatrix D mit positiven Elementen hinzu, um B = W W + D zu erhaltenk<nWk×nWWDB=WW+Dvoller Rang. Die resultierende Kovarianzmatrix kann normalisiert werden, um eine Korrelationsmatrix zu werden (wie in meiner Frage beschrieben). Das ist sehr einfach und macht den Trick. Hier sind einige beispielhafte Korrelationsmatrizen für :k=100,50,20,10,5,1

zufällige Korrelationsmatrizen aus zufälligen Faktoren

Der einzige Nachteil ist, dass die resultierende Matrix große Eigenwerte und dann einen plötzlichen Abfall aufweist, im Gegensatz zu einem guten Zerfall, der oben mit der Vine-Methode gezeigt wurde. Hier sind die entsprechenden Spektren:k

Eigenspektren dieser Matrizen

Hier ist der Code:

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)));
Amöbe sagt Reinstate Monica
quelle
+1. Hier ist jedoch nur eine Erinnerung an Ihren letzten Abschnitt über die "Faktormethode". Genau richtig ist es, wenn Spalten von Worthogonal sind (dh Kosinus zwischen ihnen ist 0). Es reicht Wnatürlich nicht aus, einfach einen Zufallsgenerator zu erstellen. Wenn sie nicht orthogonal sind - dh Faktoren sind schräg (Anruf dann Wals W_) - Faktor Satz nicht WW'aber W_CW_'mit C„Korrelationen“ (Cosinus) zwischen Faktoren zu sein. Nun, C=Q'Qmit Qder nichtorthogonalen Rotationsmatrix der Rotation W_=inv(Q)'W(und so weiter W=W_Q'). Generiere some Q- eine Matrix mit der Spalte ss = 1 und der Matrix ss = Größe der Matrix.
TTNPHNS
... Tippfehler: W_=inv(Q)'Wnatürlich nicht W_= W inv(Q)'.
ttnphns
WWW+DW
1
Übersetzen in 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
1
@Mihai, guter Punkt und deine Vorschläge sind wahrscheinlich die einfachsten. Sie könnten auchS <- matrix(nearPD(S, corr = TRUE, keepDiag = TRUE)$mat@x,ncol(S),ncol(S))
Scott Worland
7

ein

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()

Die etwas gleichmäßige Verteilung Die Ergebnisse von imshow

psarka
quelle
crsk[-ein,ein]X
Ja, du hast vollkommen recht! (Oh Junge, das war ja doof: D). Ich habe den Zufallsteil in randn (1) * a geändert und jetzt ist es viel besser.
Psarka
k
einn
Ein Nachteil dieser Methode ist, dass die resultierende Korrelationsmatrix einen großen Eigenwert aufweist, die übrigen jedoch nahezu einheitlich sind. Diese Prozedur liefert also keine "allgemeine" Korrelationsmatrix ... Nicht, dass ich sie in meiner Frage angegeben hätte. @Ssdecontrol erwähnte jedoch in den obigen Kommentaren, dass es anscheinend Möglichkeiten gibt, aus allen Korrelationsmatrizen eine Stichprobe zu erstellen. Das sieht interessant, aber viel komplizierter aus.
Amöbe sagt Reinstate Monica
6

Hmm, nachdem ich ein Beispiel in meiner MatMate-Sprache gemacht habe, sehe ich, dass es bereits eine Python-Antwort gibt, was vielleicht vorzuziehen ist, weil Python weit verbreitet ist. Aber weil Sie noch Fragen hatten, zeige ich Ihnen meinen Ansatz mit der Matmate-Matrix-Sprache, vielleicht ist er selbstkommentierender.

Methode 1
(mit 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

Das Problem hier könnte sein, dass wir Blöcke von Submatrizen definieren, die hohe Korrelationen mit geringer Korrelation aufweisen, und dies nicht programmatisch, sondern durch die konstanten Verkettungsausdrücke. Vielleicht könnte dieser Ansatz eleganter in Python modelliert werden.


Methode 2 (a)
Danach gibt es einen völlig anderen Ansatz, bei dem wir die mögliche verbleibende Kovarianz in zufälligen Mengen von 100 Prozent in eine Faktor-Ladungsmatrixfüllen. Dies geschieht in 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~

und die erzeugte Korrelationsmatrix ist

     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

Möglicherweise erzeugt dies eine Korrelationsmatrix mit dominanten Hauptkomponenten aufgrund der kumulativen Erzeugungsregel für die Faktorladungsmatrix. Es könnte auch besser sein, eine positive Bestimmtheit sicherzustellen, indem der letzte Teil der Varianz zu einem einzigartigen Faktor gemacht wird. Ich habe es im Programm belassen, um den Fokus auf das allgemeine Prinzip zu richten.

Eine 100x100-Korrelationsmatrix hatte die folgenden Korrelationshäufigkeiten (auf 1 Dezimalstelle gerundet)

    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

[aktualisieren]. Hmm, die 100x100 Matrix ist schlecht konditioniert; Pari / GP kann die Eigenwerte mit der Polroots (charpoly ()) - Funktion auch mit 200-stelliger Genauigkeit nicht korrekt ermitteln. Ich habe eine Jacobi-Rotation in pca-Form auf der Ladematrix L durchgeführt und dabei meist extrem kleine Eigenwerte gefunden, die ich in Logarithmen zur Basis 10 ausgedruckt habe (die ungefähr die Position des Dezimalpunkts angeben). Lesen Sie von links nach rechts und dann zeilenweise:

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

[Update 2]
Methode 2 (b)
Eine Verbesserung könnte darin bestehen, die artikelspezifische Varianz auf ein nicht marginales Maß zu erhöhen und auf eine einigermaßen geringere Anzahl gemeinsamer Faktoren zu reduzieren (zum Beispiel ganzzahlige Quadratwurzel der Artikelnummer):

{  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...

Die Struktur des Ergebnisses

in Bezug auf die Verteilung der Korrelationen:Bild

bleibt ähnlich (auch die böse Nicht-Zerlegbarkeit durch PariGP), aber die Eigenwerte, wenn sie durch jacobi-Rotation der Ladungsmatrix gefunden werden, haben jetzt eine bessere Struktur, für ein neu berechnetes Beispiel habe ich die Eigenwerte als erhalten

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
Gottfried Helms
quelle
Vielen Dank! Sehr interessant, aber ich brauche etwas Zeit, um es zu verdauen ...
Amöbe sagt Reinstate Monica
Ich muss Ihre Antwort noch sorgfältig durchgehen, aber in der Zwischenzeit habe ich einen Artikel über das Abtasten von Zufallskorrelationsmatrizen gelesen, und eine der Methoden von dort kann verwendet werden, um genau das zu tun, was ich brauche. Ich habe hier eine Antwort gepostet, Sie könnten interessiert sein, einen Blick darauf zu werfen! Es verweist auf eine viel detailliertere Antwort, die ich in einem anderen Thread geschrieben habe.
Amöbe sagt Reinstate Monica
@amoeba: glücklich, dass du etwas gefunden hast, das gut für dich funktioniert! Es ist eine interessante Frage, auf die ich später noch einmal zurückkommen werde, um vielleicht die MatMate-Prozeduren zu verbessern / anzupassen (und sie zu Unterprogrammen zu machen), je nachdem, an welchem ​​Papier Sie gearbeitet haben.
Gottfried Helms
2

Interessante Frage (wie immer!). Wie wäre es, eine Reihe von Beispielmatrizen zu finden, die die gewünschten Eigenschaften aufweisen, und dann konvexe Kombinationen davon zu verwenden, da wennEIN und B Sind positiv definitiv, dann ist es so λEIN+(1-λ)B. Als Bonus ist keine Neuskalierung der Diagonalen aufgrund der Konvexität der Operation erforderlich. Durch Einstellen derλUm konzentrierter gegen 0 und 1 zu sein als gleichmäßig verteilt, können Sie die Proben an den Rändern des Polytops oder im Inneren konzentrieren. (Sie können eine Beta / Dirichlet-Verteilung verwenden, um die Konzentration im Vergleich zur Gleichmäßigkeit zu steuern.)

Zum Beispiel könnten Sie lassen EIN komponentensymmetrisch sein und Bsein toeplitz. Natürlich können Sie jederzeit eine weitere Klasse hinzufügenC, und nehme λEINEIN+λBB+λCC so dass λ=1 und λ0, und so weiter.

Andrew M
quelle
Danke für den Vorschlag, Andrew, aber natürlich wäre es besser, eine unvoreingenommene Methode zu haben, die nicht mit vordefinierten Methoden beginnen muss EIN und B... In den Kommentaren zu meinem ursprünglichen Frage @ssdecontrol auf ein Papier bezeichnet Algorithmen zur Probe Korrelationsmatrizen einheitlich zu beschreiben (in einem bestimmten präzisen Sinne) oder in Richtung Einheitsmatrix voreingenommen, aber ich kann nicht einen Weg finden , um sie noch voreingenommen Probe weg von Identität ... Ich habe hier auch ein paar alte Threads gefunden, die fast die gleiche Frage gestellt haben. Vielleicht interessiert Sie das, sehen Sie sich mein Update an.
Amöbe sagt Reinstate Monica
Ah, aber mit einem solchen Algorithmus und einer geeigneten Diversität der "Eckpunkte" (dh Matrizen), die Ihr Polytop aus positiv-definitiven Korrelationsmatrizen definieren, können Sie die Ablehnungsabtastung verwenden, um die Verteilung der Eigenwerte und die Gleichmäßigkeit der Einträge zu ermitteln. usw., die Sie wünschen. Mir ist jedoch nicht klar, was eine gute Basis wäre. Klingt nach einer Frage für jemanden, der in jüngerer Zeit abstrakte Algebra studiert hat als ich.
Andrew M
Hallo zusammen, ich habe einen Artikel über das Abtasten von Zufallskorrelationsmatrizen gelesen und eine der Methoden von dort kann verwendet werden, um genau das zu tun, was ich brauche. Ich habe hier eine Antwort gepostet, Sie könnten interessiert sein, einen Blick darauf zu werfen! Es verweist auf eine viel detailliertere Antwort, die ich in einem anderen Thread geschrieben habe.
Amöbe sagt Reinstate Monica
2

R hat ein Paket (clusterGeneration) , das die Methode implementiert in:

Beispiel:

> (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

Leider scheint es nicht möglich zu sein, Korrelationen zu simulieren, die einer gleichmäßigen Verteilung folgen. Es scheint stärkere Korrelationen herzustellen, wenn alphadauf sehr kleine Werte gesetzt wird, aber selbst bei würde 1/100000000000000der Bereich der Korrelationen nur bis zu ungefähr 1,40 reichen.

Trotzdem hoffe ich, dass dies für jemanden von Nutzen sein könnte.

Deleet
quelle