Cara menggambar plot scree dengan python? [Tutup]

11

Saya menggunakan dekomposisi vektor singular pada sebuah matriks dan mendapatkan matriks U, S dan Vt. Pada titik ini, saya mencoba untuk memilih ambang batas untuk mempertahankan jumlah dimensi. Saya disarankan untuk melihat plot scree tetapi bertanya-tanya bagaimana cara merencanakannya dalam numpy. Saat ini, saya melakukan hal berikut menggunakan perpustakaan numpy dan scipy dengan python:

U, S, Vt = svd(A)

Ada saran?

Legenda
sumber
1
ambil diagonal S, jika belum menjadi diagonal, kuadratkan, urutkan dalam urutan menurun, ambil jumlah kumulatif, bagi dengan nilai terakhir, lalu plot.
shabbychef
@shabbychef: Maksud Anda, ambil jumlah kumulatif dan bagi dengan semua nilai, kan?
Legenda
Iya. Di matlab, itu akan menjadi[U,S,V] = svd(X);S = cumsum(sort(diag(S).^2,1,'descend'));S = S ./ S(end);plot(S);
shabbychef

Jawaban:

13

Berikut ini adalah contoh yang dapat disisipkan ke prompt IPython dan menghasilkan gambar seperti di bawah ini (menggunakan data acak):

import numpy as np
import matplotlib
import matplotlib.pyplot as plt

#Make a random array and then make it positive-definite
num_vars = 6
num_obs = 9
A = np.random.randn(num_obs, num_vars)
A = np.asmatrix(A.T) * np.asmatrix(A)
U, S, V = np.linalg.svd(A) 
eigvals = S**2 / np.sum(S**2)  # NOTE (@amoeba): These are not PCA eigenvalues. 
                               # This question is about SVD.

fig = plt.figure(figsize=(8,5))
sing_vals = np.arange(num_vars) + 1
plt.plot(sing_vals, eigvals, 'ro-', linewidth=2)
plt.title('Scree Plot')
plt.xlabel('Principal Component')
plt.ylabel('Eigenvalue')
#I don't like the default legend so I typically make mine like below, e.g.
#with smaller fonts and a bit transparent so I do not cover up data, and make
#it moveable by the viewer in case upper-right is a bad place for it 
leg = plt.legend(['Eigenvalues from SVD'], loc='best', borderpad=0.3, 
                 shadow=False, prop=matplotlib.font_manager.FontProperties(size='small'),
                 markerscale=0.4)
leg.get_frame().set_alpha(0.4)
leg.draggable(state=True)
plt.show()

masukkan deskripsi gambar di sini

Josh Hemann
sumber
Hermann: +1 Terima kasih atas waktu Anda! Saya tahu ini sudah lama tetapi bagaimanapun ini benar-benar baik untuk dimiliki :)
Legenda
apa num_vars? sepertinya tidak didefinisikan dalam skrip Anda.
TheChymera
@TheChymera - Terima kasih telah menangkap ini, saya telah memperbarui tanggapan saya.
Josh Hemann
@Josh Hemann ya, saya juga menemukan ini dalam waktu yang bersamaan - tapi saya pikir mungkin lebih baik untuk menghitungnya dari bentuk A
TheChymera
1
@ JoshHemann dapatkah Anda menjelaskan ini: eigvals = S ** 2 / np.cumsum (S) [- 1] ?? Saya telah melihat berdasarkan pada beberapa makalah eigvals = S ** 2 / (n-1) di mana n adalah jumlah fitur
makis