Come disegnare campioni dalla distribuzione gaussiana multivariata parametrizzata con precisione in MATLAB

Question

Mi chiedo come disegnare campioni in MATLAB, dove ho matrice di precisione e intendo come argomento di input.

So che mvnrnd è un modo tipico per farlo, ma richiede la matrice di covarianza (cioè inverso di precisione)) come argomento.

Ho solo avere matrice di precisione, ed a causa della questione di calcolo, non riesco a invertire la mia matrice di precisione, dal momento che ci vorrà troppo tempo (la mia dimensione è di circa 2000 * 2000)

Answer 1

Buona domanda. Si noti che è possibile generare campioni da una distribuzione normale multivariante utilizzando campioni della distribuzione normale standard mediante la procedura descritta in the relevant Wikipedia article.

Fondamentalmente, questo si riduce a valutare A*z + mu dove z è un vettore di variabili casuali indipendenti campione dalla distribuzione normale standard, mu è un vettore di mezzo, e A*A' = Sigma è la matrice di covarianza. Poiché si ha l'inverso di quest'ultima quantità, ovvero inv(Sigma), è possibile eseguire probabilmente una scomposizione di Cholesky (vedere chol) per determinare l'inverso di A. È quindi necessario valutare A * z. Se si conosce solo inv(A), ciò può ancora essere eseguito senza eseguire una matrice inversa risolvendo invece un sistema lineare (ad esempio tramite l'operatore backslash).

La scomposizione di Cholesky potrebbe essere ancora problematica per voi, ma spero che questo aiuti.

Answer 2

Se volete assaggiare da N (μ, Q ^-1) e solo Q è disponibile, si può prendere la fattorizzazione di Cholesky di Q, L, in modo tale che LL ^T = Q. Quindi prendere l'inverso di L ^T, L ^-T e campionare Z da una normale distribuzione normale N (0, I).

Considerando che L ^-T è una matrice triangolare superiore dxd e Z è un vettore colonna d-dimensionale, μ + L ^-T Z sarà distribuito come N (μ, Q ^-1) .

Se si desidera evitare l'inverso di L, è invece possibile risolvere il sistema triangolare di equazioni L ^T v = Z mediante sostituzione posteriore. μ + v sarà quindi distribuito come N (μ, Q ^-1).

Alcuni codice MATLAB illustrativi:

% fa una matrice 2x2 di covarianza e un vettore medio

covm = [3 0.4 * (sqrt (3 * 7)); 0.4 * (sqrt (3 * 7)) 7];

mu = [100; 2];

% Ottenere la matrice precisione

Q = inv (covm);

% prendere la decomposizione di Cholesky di Q (Chol in MATLAB restituisce già il fattore triangolare superiore)

L = Chol (Q);

% Pari 2000 campioni provenienti da una distribuzione normale bivariata normale

Z = normrnd (0,1, [2, 2000]);

% risolvere il sistema e aggiungere la media

X = repmat (mu, 1, 2000) + L \ Z;

% controllare il risultato

medio (X ')

var (X')

corrcoef (X ')

% confrontare il campionamento dalla matrice di covarianza

Y = mvnrnd (mu, covm, 2000) ';

medio (Y ')

var (Y')

corrcoef (Y ')

dispersione (X (1, ​​:), X (2, :),' b')

presa su

scatter (Y (1, :), Y (2, :), 'r')

per una maggiore efficienza, credo che è possibile cercare som Il pacchetto che risolve in modo efficiente i sistemi triangolari.