Как нарисовать осыпь в Python? [закрыто]

11

Я использую сингулярное векторное разложение на матрице и получаю матрицы U, S и Vt. На данный момент я пытаюсь выбрать порог для количества измерений, чтобы сохранить. Мне предложили взглянуть на сюжетный осколок, но мне интересно, как сделать так, чтобы он наносился на ноль. В настоящее время я делаю следующее, используя библиотеки numpy и scipy в python:

U, S, Vt = svd(A)

Какие-либо предложения?

data-visualization python svd легенда
источник
1
возьмите диагональ S, если она еще не диагональ, возведите ее в квадрат, отсортируйте в порядке убывания, возьмите накопленную сумму, разделите на последнее значение, затем построите ее.
Шаббычеф
@shabbychef: Вы имеете в виду, возьмите совокупную сумму и поделите на сумму всех значений, верно?
Легенда
да. В Matlab это было бы[U,S,V] = svd(X);S = cumsum(sort(diag(S).^2,1,'descend'));S = S ./ S(end);plot(S);
шаббучеф

Ответы:

13

Вот пример, который можно вставить в приглашение IPython и создать изображение, как показано ниже (оно использует случайные данные):

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

введите описание изображения здесь

Джош Хеманн
источник
Германн: +1 Спасибо за потраченное время! Я знаю, что это было давно, но, тем не менее, это действительно хорошо иметь :)
Legend
что такое num_vars? кажется, это не определено в вашем сценарии.
TheChemera
@TheChymera - Спасибо, что поймали это, я обновил свой ответ.
Джош Хеманн
@ Джош Хеманн: Да, я также понял это в то же время - но я думаю, что было бы лучше вычислить его по форме A
TheChymera
1
@JoshHemann Вы можете объяснить это: eigvals = S ** 2 / np.cumsum (S) [- 1] ?? На основании некоторых работ я видел eigvals = S ** 2 / (n-1), где n - количество функций
makis