SVD, pseudoinverse en PCA: SVD, pseudoinverse en PCA in MATLAB
PCA
Laten we in een concreet voorbeeld de principale componenten analyse in MATLAB uitvoeren.
We beginnen met een dataset te creëren voor twee grootheden en de gegevens te visualiseren.
>> rng(0); % controleer random number generatie
>> data(:,1) = randn(25,1);
>> data(:,2) = 2 + 1.5*data(:,1) + 0.25*randn(size(data(:,1))); % data op een lijn + ruis
>> data = sortrows(data,1); % sorteer eerste kolom
>> figure
axes('LineWidth',0.6,...
'FontName','Helvetica',...
'FontSize',8,...
'XAxisLocation','Origin',...
'YAxisLocation','Origin');
line(data(:,1),data(:,2),...
'LineStyle','None',...
'Marker','o');
axis equal
Om de analyse goed te kunnen uitvoeren moeten we de gegevens eerst bewerken zodanig dat de gemiddelde waarden van de twee grootheden gelijk aan nul zijn.
>> data(:,1) = data(:,1)-mean(data(:,1));
>> data(:,2) = data(:,2)-mean(data(:,2));
We kunnen voor de datamatrix \(X\) met 2 grootheden (2 kolommen) en \(N=25\) gegevens (25 rijen) dan de covariantiematrix \(C\) berekenen als \[C = \frac{1}{N-1}X^{\top}\!X\] MATLAB heeft hiervoor de functie cov
aan boord.
>> C = cov(data)
C =
1.9179 2.7857
2.7857 4.1102
We willen nu nieuwe grootheden definiëren wier waarden ongecorreleerd zijn. Hiervoor diagonaliseren we de covariantiematrix door het eigenwaardenprobleem voor de matrix \(C\) numeriek op te lossen.
>> [V,D] = eig(C)
V =
-0.8265 0.5630
0.5630 0.8265
D =
0.0205 0
0 6.0076
De eigenvectoren definiëren nieuwe grootheden (als lineaire combinaties van de oude grootheden) die ongecorreleerd zijn omdat correlatie niets ander is dan inproduct gelijk aan nul voor de kolomvectoren. In onderstaande figuur staan de twee vectoren afgebeeld als stukje van een nieuw coordinatensysteem.
>> figure
axes('LineWidth',0.6,...
'FontName','Helvetica',...
'FontSize',8,...
'XAxisLocation','Origin',...
'YAxisLocation','Origin');
line(data(:,1),data(:,2),...
'LineStyle','None',...
'Marker','o');
line([0 V(1,1)],[0 V(2,1)],...
'Color',[0.8 0.5 0.3],...
'LineWidth',0.75);
line([0 V(1,2)],[0 V(2,2)],...
'Color',[0.8 0.5 0.3],...
'LineWidth',0.75);
axis equal
We kunnen de gegevens transformeren voor het nieuwe coördinatensysteem en visualiseren. Omdat de eigenwaarden in opklimmende volgorde staan moeten we tussendoor kolommen verwisselen.
>> newdata = data * V';
>> newdata = fliplr(newdata)
newdata =
-5.2475 -0.0866
-3.2490 0.1749
-3.2866 0.0982
-3.3236 -0.0482
-2.1315 -0.1726
-1.5034 -0.0214
-1.3865 -0.0396
-0.9434 0.1884
-0.7009 -0.1084
-0.7160 -0.1476
-0.0023 0.1615
0.0702 0.1519
0.0915 0.0045
-0.1780 -0.2314
0.3098 0.0978
0.0389 -0.0966
0.5021 0.0535
1.4919 0.0659
1.2632 -0.0997
1.6386 0.0683
2.0519 0.1797
2.3433 0.1318
3.2713 -0.3681
4.4254 0.0968
5.1706 -0.0531
>> figure
axes('LineWidth',0.6,...
'FontName','Helvetica',...
'FontSize',8,...
'XAxisLocation','Origin',...
'YAxisLocation','Origin')
line(newdata(:,1),newdata(:,2),...
'LineStyle','None',...
'Marker','o');
axis equal
Met de functie var
kunnen we de variatie in de grootheden bepalen. De ene grootheid betreft dan \(99.96\%\) van de variatie van de data, en de andere grootheid slechts voor \(0.34\%\).
>> var(newdata)
ans =
6.0076 0.0205
>> var(newdata)/sum(var(newdata))
ans =
0.9966 0.0034
\(\phantom{x}\)
We kunnen natuurlijk ook de ingebouwde functie pca
in MATLAB gebruiken.
>> [coeff,newdata,latent,tsd,variance] = pca(data);
>> coeff % coëfficiënten van de principale componenten (factorladingen)
coeff =
0.5630 0.8265
0.8265 -0.5630
>> newdata % getransformeerde data
newdata =
-5.2475 0.0866
-3.2490 -0.1749
-3.2866 -0.0982
-3.3236 0.0482
-2.1315 0.1726
-1.5034 0.0214
-1.3865 0.0396
-0.9434 -0.1884
-0.7009 0.1084
-0.7160 0.1476
-0.0023 -0.1615
0.0702 -0.1519
0.0915 -0.0045
-0.1780 0.2314
0.3098 -0.0978
0.0389 0.0966
0.5021 -0.0535
1.4919 -0.0659
1.2632 0.0997
1.6386 -0.0683
2.0519 -0.1797
2.3433 -0.1318
3.2713 0.3681
4.4254 -0.0968
5.1706 0.0531
>> latent % eigenwaarden van de covariantiematrix bij de dataset
latent =
6.0076
0.0205
>> variance % variatie in data langs de principale assen
variance =
99.6601
0.3399
>> corrcoef(newdata)
ans =
1.0000 0.0000
0.0000 1.0000