Cara membuat matriks kovarians yang sewenang-wenang

21

Misalnya, dalam R, MASS::mvrnorm()fungsi ini berguna untuk menghasilkan data untuk menunjukkan berbagai hal dalam statistik. Dibutuhkan Sigmaargumen wajib yang merupakan matriks simetris yang menentukan matriks kovarians dari variabel. Bagaimana saya membuat simetris n×n matriks dengan entri sewenang-wenang?

rsl
sumber
3
Saya pikir pertanyaan ini akan mendapat manfaat dari diedit untuk fokus pada "bagaimana saya bisa membuat matriks kovarian sewenang-wenang" dan kurang pada aspek pengkodean. Tentu saja ada masalah statistik yang mendasari topik di sini, seperti yang ditunjukkan oleh jawabannya.
Silverfish

Jawaban:

22

Buat matriks A dengan nilai arbitrern×nSEBUAH

dan kemudian gunakan sebagai matriks kovarians Anda. Σ=SEBUAHTSEBUAH

Sebagai contoh

n <- 4  
A <- matrix(runif(n^2)*2-1, ncol=n) 
Sigma <- t(A) %*% A
Henry
sumber
Demikian juga Sigma <- A + t(A),.
rsl
6
@MoazzemHossen: Saran Anda akan menghasilkan matriks simetris, tetapi mungkin tidak selalu berupa semidefinit positif (mis. Saran Anda dapat menghasilkan matriks dengan nilai eigen negatif) dan karenanya mungkin tidak cocok sebagai matriks kovarians
Henry
Ya, saya perhatikan bahwa R mengembalikan kesalahan jika cara yang saya sarankan menghasilkan matriks yang tidak sesuai.
rsl
4
Perhatikan bahwa jika Anda lebih suka matriks korelasi untuk interpretabilitas yang lebih baik, ada fungsi ? Cov2cor , yang dapat diterapkan selanjutnya.
gung - Reinstate Monica
1
@ B11b: Anda perlu matriks kovarians Anda untuk menjadi semi-pasti positif. Itu akan membatasi nilai kovarian, tidak sepenuhnya jelas ketika n>2
Henry
24

Saya suka memiliki kontrol atas objek yang saya buat, bahkan ketika mereka mungkin sewenang-wenang.

Pertimbangkan, kemudian, bahwa semua kemungkinan kovarians matriks Σ dapat dinyatakan dalam bentukn×nΣ

Σ=P Diagonal(σ1,σ2,,σn) P

di mana adalah matriks ortogonal dan σ 1σ 2σ n0 .Pσ1σ2σn0

Secara geometris ini menggambarkan struktur kovarians dengan berbagai komponen utama ukuran . Komponen ini menunjukkan dalam arah dari barisan P . Lihat gambar di Membuat analisis komponen utama, vektor eigen & nilai eigen untuk contoh dengan n = 3 . Pengaturan σ i akan mengatur besaran kovarian dan ukuran relatifnya, sehingga menentukan bentuk ellipsoidal yang diinginkan. Baris P mengorientasikan sumbu bentuk yang Anda inginkan.σiPn=3σiP

Salah satu manfaat aljabar dan komputasi dari pendekatan ini adalah ketika , Σ mudah terbalik (yang merupakan operasi umum pada matriks kovarians):σn>0Σ

Σ1=P Diagonal(1/σ1,1/σ2,,1/σn) P.

Tidak peduli dengan arah, tetapi hanya tentang rentang ukuran ? Tidak apa-apa: Anda dapat dengan mudah menghasilkan matriks ortogonal acak. Hanya membungkus n 2 nilai standar iid Normal ke dalam matriks persegi dan kemudian orthogonalize itu. Ini hampir pasti akan berhasil (asalkan n tidak besar). Dekomposisi QR akan melakukan itu, seperti dalam kode iniσin2n

n <- 5
p <- qr.Q(qr(matrix(rnorm(n^2), n)))

Ini berfungsi karena distribusi multinormal -variat yang dihasilkan adalah "elips": ia tidak berubah dalam semua rotasi dan refleksi (melalui titik asal). Jadi, semua matriks ortogonal dihasilkan secara seragam, seperti yang diperdebatkan di Bagaimana menghasilkan titik-titik yang terdistribusi secara merata pada permukaan bola unit 3-d? .n

Cara cepat untuk mendapatkan dari P dan σ i , setelah Anda menentukan atau membuat mereka, menggunakan dan mengeksploitasi penggunaan kembali array dalam operasi aritmatika, seperti dalam contoh ini dengan σ = ( σ 1 , , σ 5 ) = ( 5 , 4 , 3 , 2 , 1 ) :ΣPσicrossprodRσ=(σ1,,σ5)=(5,4,3,2,1)

Sigma <- crossprod(p, p*(5:1))

Sebagai tanda centang, dekomposisi Nilai Singular harus mengembalikan dan P . Anda dapat memeriksanya dengan perintahσP

svd(Sigma)

Sigmaσ

Tau <- crossprod(p, p/(5:1))

zapsmall(Sigma %*% Tau)n×nσi01/σiσi

whuber
sumber
P
1
Mungkin perlu disebutkan bahwa nilai-nilai tunggal di svd(Sigma)akan ditata ulang - yang membingungkan saya sebentar.
FrankD
1

Anda dapat mensimulasikan matriks definitif positif acak dari distribusi Wishart menggunakan fungsi "rWishart" dari paket "stats".

n <- 4
rWishart(1,n,diag(n))
Carlos Llosa
sumber
1

Ada paket khusus untuk itu, clusterGeneration(ditulis antara lain oleh Harry Joe, nama besar di bidang itu).

Ada dua fungsi utama:

  • genPositiveDefMat menghasilkan matriks kovarians, 4 metode berbeda
  • rcorrmatrix : menghasilkan matriks korelasi

Contoh cepat:

library(clusterGeneration)
#> Loading required package: MASS
genPositiveDefMat("unifcorrmat",dim=3)
#> $egvalues
#> [1] 15.408962  5.673916  1.228842
#> 
#> $Sigma
#>          [,1]     [,2]     [,3]
#> [1,] 6.714871 1.643449 6.530493
#> [2,] 1.643449 6.568033 2.312455
#> [3,] 6.530493 2.312455 9.028815
genPositiveDefMat("eigen",dim=3)
#> $egvalues
#> [1] 8.409136 4.076442 2.256715
#> 
#> $Sigma
#>            [,1]       [,2]      [,3]
#> [1,]  2.3217300 -0.1467812 0.5220522
#> [2,] -0.1467812  4.1126757 0.5049819
#> [3,]  0.5220522  0.5049819 8.3078880

Dibuat pada 2019-10-27 oleh paket reprex (v0.3.0)

Akhirnya, perhatikan bahwa pendekatan alternatif adalah melakukan percobaan pertama dari awal, kemudian gunakan Matrix::nearPD()untuk membuat matriks Anda pasti-positif.

Matifou
sumber