2013-08-18 14 views
7

Próbuję wykonać prostą analizę głównych składników za pomocą matplotlib.mlab.PCA, ale z atrybutami klasy nie mogę uzyskać czystego rozwiązania mojego problemu. Oto przykład:Podstawowy przykład dla PCA z matplotlib

Get niektóre dane manekina w 2D i rozpocząć PCA:

from matplotlib.mlab import PCA 
import numpy as np 

N  = 1000 
xTrue = np.linspace(0,1000,N) 
yTrue = 3*xTrue 

xData = xTrue + np.random.normal(0, 100, N) 
yData = yTrue + np.random.normal(0, 100, N) 
xData = np.reshape(xData, (N, 1)) 
yData = np.reshape(yData, (N, 1)) 
data = np.hstack((xData, yData)) 
test2PCA = PCA(data) 

Teraz, po prostu chcę, aby uzyskać podstawowe komponenty jako wektory w moich oryginalnych współrzędnych i wykreślić je jako strzałki na przetwarzanie moich danych.

Jaki jest szybki i czysty sposób, aby się tam dostać?

Dzięki Tyrax

Odpowiedz

22

Nie sądzę klasa mlab.PCA jest odpowiednia do tego, co chcesz zrobić. W szczególności, klasa PCA przeskalowanie danych przed znalezienie wektorów własnych:

a = self.center(a) 
U, s, Vh = np.linalg.svd(a, full_matrices=False) 

Sposób center dzieli przez sigma:

def center(self, x): 
    'center the data using the mean and sigma from training set a' 
    return (x - self.mu)/self.sigma 

Powoduje to wektory własne, pca.Wt, tak:

[[-0.70710678 -0.70710678] 
[-0.70710678 0.70710678]] 

Są prostopadłe, ale nie odnoszą się bezpośrednio do głównych osi oryginalnych danych. Są to główne osie w odniesieniu do masowanych danych.

Być może łatwiej będzie kodować co chcesz bezpośrednio (bez użycia klasy mlab.PCA):

import numpy as np 
import matplotlib.pyplot as plt 

N = 1000 
xTrue = np.linspace(0, 1000, N) 
yTrue = 3 * xTrue 
xData = xTrue + np.random.normal(0, 100, N) 
yData = yTrue + np.random.normal(0, 100, N) 
xData = np.reshape(xData, (N, 1)) 
yData = np.reshape(yData, (N, 1)) 
data = np.hstack((xData, yData)) 

mu = data.mean(axis=0) 
data = data - mu 
# data = (data - mu)/data.std(axis=0) # Uncommenting this reproduces mlab.PCA results 
eigenvectors, eigenvalues, V = np.linalg.svd(data.T, full_matrices=False) 
projected_data = np.dot(data, eigenvectors) 
sigma = projected_data.std(axis=0).mean() 
print(eigenvectors) 

fig, ax = plt.subplots() 
ax.scatter(xData, yData) 
for axis in eigenvectors: 
    start, end = mu, mu + sigma * axis 
    ax.annotate(
     '', xy=end, xycoords='data', 
     xytext=start, textcoords='data', 
     arrowprops=dict(facecolor='red', width=2.0)) 
ax.set_aspect('equal') 
plt.show() 

enter image description here

+0

wielkie dzięki. Tego właśnie szukałem. – Tyrax

+0

jakie znaczenie ma stała 1.618? Skąd to pochodzi ? – joaquin

+1

@joaquin: Jest to w przybliżeniu [złoty współczynnik] (http://en.wikipedia.org/wiki/Golden_ratio). Możesz oczywiście wybrać dowolną stałą, którą lubisz, ale [często wygląda dobrze] (http://en.wikipedia.org/wiki/Golden_ratio#Painting). – unutbu

Powiązane problemy