pca
生データの主成分分析
構文
説明
coeff = pca(X,Name,Value)Name,Value のペア引数により指定され、特殊なデータ型の計算と処理を行う追加オプションを使用すると、前の構文に任意の出力引数を返します。
たとえば、pca から返される主成分の数または使用する SVD 以外のアルゴリズムを指定することができます。
例
標本データ セットを読み込みます。
load hald成分データには、4 つの変数に対しての 13 の観測値があります。
成分データの主成分を求めます。
coeff = pca(ingredients)
coeff = 4×4
-0.0678 -0.6460 0.5673 0.5062
-0.6785 -0.0200 -0.5440 0.4933
0.0290 0.7553 0.4036 0.5156
0.7309 -0.1085 -0.4684 0.4844
coeff の行には 4 つの変数の係数が含まれ、その列は 4 つの主成分に一致します。
データ セットに欠損値がある場合に主成分係数を求めます。
標本データ セットを読み込みます。
load imports-85データ行列 X には、3 列目から 15 列目に 13 個の連続変数が含まれています (wheel-base、length、width、height、curb-weight、engine-size、bore、stroke、compression-ratio、horsepower、peak-rpm、city-mpg、highway-mpg)。bore と stroke の各変数の 56 ～ 59 行目の 4 つの値と、horsepower と peak-rpm の各変数の 131 行目と 132 行目の 2 つの値が欠けています。
主成分分析を実行します。
coeff = pca(X(:,3:15));
既定では pca が 'Rows','complete' の名前と値のペア引数によって指定されたアクションを実行します。このオプションは計算前に NaN の値をもつ観測を削除します。score および tsquared の対応する位置、つまり 56 ～ 59 行目、131 行目および 132 行目には NaN の行が再挿入されます。
'pairwise' を使用して主成分分析を実行します。
coeff = pca(X(:,3:15),'Rows','pairwise');
ここでは pca が、X の i 列または j 列に NaN 値のない行を使用して、共分散行列の (i, j) 要素を計算します。このとき、結果の共分散行列は正定値でない場合があります。このオプションを使用できるのは、pca が固有値分解アルゴリズムを使用する場合のみです。この例のようにアルゴリズムを指定しない場合、アルゴリズムは pca によって 'eig' に設定されます。'pairwise' オプションとともにアルゴリズムとして 'svd' を指定すると、pca は警告メッセージを返し、アルゴリズムを 'eig' に設定してから処理を続行します。
'Rows','all' の名前と値のペア引数を使用した場合、このオプションではデータ セットに欠損値がないと仮定されるため、pca は終了します。
coeff = pca(X(:,3:15),'Rows','all');
Error using pca (line 180) Raw data contains NaN missing value while 'Rows' option is set to 'all'. Consider using 'complete' or pairwise' option instead.
主成分分析の実行中に変数の逆分散を重みとして使用します。
標本データ セットを読み込みます。
load hald成分の分散の逆数を変数の重みとして使用して、主成分分析を実行します。
[wcoeff,~,latent,~,explained] = pca(ingredients,'VariableWeights','variance')
wcoeff = 4×4
-2.7998 2.9940 -3.9736 1.4180
-8.7743 -6.4411 4.8927 9.9863
2.5240 -3.8749 -4.0845 1.7196
9.1714 7.5529 3.2710 11.3273
latent = 4×1
2.2357
1.5761
0.1866
0.0016
explained = 4×1
55.8926
39.4017
4.6652
0.0406
係数行列 wcoeff は正規直交ではない点に注意してください。
正規直交係数行列を計算します。
coefforth = diag(std(ingredients))\wcoeff
coefforth = 4×4
-0.4760 0.5090 -0.6755 0.2411
-0.5639 -0.4139 0.3144 0.6418
0.3941 -0.6050 -0.6377 0.2685
0.5479 0.4512 0.1954 0.6767
新しい係数行列 coefforth の正規直交性を確認します。
coefforth*coefforth'
ans = 4×4
1.0000 0.0000 -0.0000 0.0000
0.0000 1.0000 0.0000 0
-0.0000 0.0000 1.0000 -0.0000
0.0000 0 -0.0000 1.0000
データに欠損値があるときに交互最小二乗 (ALS) アルゴリズムを使用して主成分を見つけます。
標本データを読み込みます。
load hald成分データには、4 つの変数に対しての 13 の観測値があります。
ALS アルゴリズムを使って主成分分析を実行し、成分係数を表示します。
[coeff,score,latent,tsquared,explained] = pca(ingredients); coeff
coeff = 4×4
-0.0678 -0.6460 0.5673 0.5062
-0.6785 -0.0200 -0.5440 0.4933
0.0290 0.7553 0.4036 0.5156
0.7309 -0.1085 -0.4684 0.4844
欠損値を無作為に導入します。
y = ingredients; rng('default'); % for reproducibility ix = random('unif',0,1,size(y))<0.30; y(ix) = NaN
y = 13×4
7 26 6 NaN
1 29 15 52
NaN NaN 8 20
11 31 NaN 47
7 52 6 33
NaN 55 NaN NaN
NaN 71 NaN 6
1 31 NaN 44
2 NaN NaN 22
21 47 4 26
NaN 40 23 34
11 66 9 NaN
10 68 8 12
⋮
データのおよそ 30% に、NaN で示される欠損値があります。
ALS アルゴリズムを使って主成分分析を実行し、成分係数を表示します。
[coeff1,score1,latent,tsquared,explained,mu1] = pca(y,... 'algorithm','als'); coeff1
coeff1 = 4×4
-0.0362 0.8215 -0.5252 0.2190
-0.6831 -0.0998 0.1828 0.6999
0.0169 0.5575 0.8215 -0.1185
0.7292 -0.0657 0.1261 0.6694
推定平均値を表示します。
mu1
mu1 = 1×4
8.9956 47.9088 9.0451 28.5515
観測されたデータを再構成します。
t = score1*coeff1' + repmat(mu1,13,1)
t = 13×4
7.0000 26.0000 6.0000 51.5250
1.0000 29.0000 15.0000 52.0000
10.7819 53.0230 8.0000 20.0000
11.0000 31.0000 13.5500 47.0000
7.0000 52.0000 6.0000 33.0000
10.4818 55.0000 7.8328 17.9362
3.0982 71.0000 11.9491 6.0000
1.0000 31.0000 -0.5161 44.0000
2.0000 53.7914 5.7710 22.0000
21.0000 47.0000 4.0000 26.0000
21.5809 40.0000 23.0000 34.0000
11.0000 66.0000 9.0000 5.7078
10.0000 68.0000 8.0000 12.0000
⋮
ALS アルゴリズムで、データの欠損値が推定されます。
結果を比較するもう 1 つの方法として、係数ベクトルの範囲内の 2 つの空間に存在する角度を見つけます。ALS を使用して、完全なデータで見つかった係数と欠損値をもつデータで見つかった係数間の角度を見つけます。
subspace(coeff,coeff1)
ans = 8.7537e-16
これは小さい値です。これは、欠損データがないときに 'Rows','complete' の名前と値のペア引数を含む pca を使用する結果と、欠損データがあるときに 'algorithm','als' の名前と値のペア引数を含む pca を使用する結果が、相互に近いことを示します。
'Rows','complete' の名前と値のペア引数を使用して主成分分析を実行し、成分係数を表示します。
[coeff2,score2,latent,tsquared,explained,mu2] = pca(y,... 'Rows','complete'); coeff2
coeff2 = 4×3
-0.2054 0.8587 0.0492
-0.6694 -0.3720 0.5510
0.1474 -0.3513 -0.5187
0.6986 -0.0298 0.6518
この例で、pca が欠損値をもつ行を削除するため、y には欠損値のない 4 つの行だけが含まれるようになります。pca は 3 つの主成分のみを返します。共分散行列は半正定値ではなく、pca からエラー メッセージが返されるため、'Rows','pairwise' オプションは使用できません。
完全なデータで見つかった係数と、リストワイズ除去 ('Rows','complete' の場合) で欠損値をもつデータで見つかった係数との間の角度を検出します。
subspace(coeff(:,1:3),coeff2)
ans = 0.3576
2 つの空間の角度はかなり大きくなります。これは 2 つの結果が異なることを示しています。
推定平均値を表示します。
mu2
mu2 = 1×4
7.8889 46.9091 9.8750 29.6000
この場合、平均値は単なる y の標本平均です。
観測されたデータを再構成します。
score2*coeff2'
ans = 13×4
NaN NaN NaN NaN
-7.5162 -18.3545 4.0968 22.0056
NaN NaN NaN NaN
NaN NaN NaN NaN
-0.5644 5.3213 -3.3432 3.6040
NaN NaN NaN NaN
NaN NaN NaN NaN
NaN NaN NaN NaN
NaN NaN NaN NaN
12.8315 -0.1076 -6.3333 -3.7758
NaN NaN NaN NaN
NaN NaN NaN NaN
1.4680 20.6342 -2.9292 -18.0043
⋮
これは、NaN 値を含む行の削除が、ALS アルゴリズムのように機能しないことを示しています。データの欠損値が多すぎる場合は、ALS を使用する方が適切です。
主成分の係数、スコアおよび分散を求めます。
標本データ セットを読み込みます。
load hald成分データには、4 つの変数に対しての 13 の観測値があります。
主成分の成分データの係数、スコアおよび分散を計算します。
[coeff,score,latent] = pca(ingredients)
coeff = 4×4
-0.0678 -0.6460 0.5673 0.5062
-0.6785 -0.0200 -0.5440 0.4933
0.0290 0.7553 0.4036 0.5156
0.7309 -0.1085 -0.4684 0.4844
score = 13×4
36.8218 -6.8709 -4.5909 0.3967
29.6073 4.6109 -2.2476 -0.3958
-12.9818 -4.2049 0.9022 -1.1261
23.7147 -6.6341 1.8547 -0.3786
-0.5532 -4.4617 -6.0874 0.1424
-10.8125 -3.6466 0.9130 -0.1350
-32.5882 8.9798 -1.6063 0.0818
22.6064 10.7259 3.2365 0.3243
-9.2626 8.9854 -0.0169 -0.5437
-3.2840 -14.1573 7.0465 0.3405
9.2200 12.3861 3.4283 0.4352
-25.5849 -2.7817 -0.3867 0.4468
-26.9032 -2.9310 -2.4455 0.4116
⋮
latent = 4×1
517.7969
67.4964
12.4054
0.2372
score の各列が 1 つの主成分に対応します。ベクトル latent には 4 つの主成分の分散が保存されます。
センタリングされた成分データを復元します。
Xcentered = score*coeff'
Xcentered = 13×4
-0.4615 -22.1538 -5.7692 30.0000
-6.4615 -19.1538 3.2308 22.0000
3.5385 7.8462 -3.7692 -10.0000
3.5385 -17.1538 -3.7692 17.0000
-0.4615 3.8462 -5.7692 3.0000
3.5385 6.8462 -2.7692 -8.0000
-4.4615 22.8462 5.2308 -24.0000
-6.4615 -17.1538 10.2308 14.0000
-5.4615 5.8462 6.2308 -8.0000
13.5385 -1.1538 -7.7692 -4.0000
-6.4615 -8.1538 11.2308 4.0000
3.5385 17.8462 -2.7692 -18.0000
2.5385 19.8462 -3.7692 -18.0000
⋮
Xcentered の新しいデータは、対応する列から列平均を差し引くことで元の成分データがセンタリングされています。
各変数の正規直交主成分係数と各観測の主成分得点を 1 つのプロットで可視化します。
biplot(coeff(:,1:2),'scores',score(:,1:2),'varlabels',{'v_1','v_2','v_3','v_4'});
このバイプロットでは、4 つの変数がすべてベクトルで表されます。ベクトルの向きと大きさは、プロットにおける 2 つの主成分に対する各変数の寄与の程度を表します。たとえば、横軸上にある最初の主成分には、3 番目の変数と 4 番目の変数に対する正の係数があります。したがって、ベクトル および はプロットの右半分に向いています。1 番目の主成分で最も大きい係数は、変数 に対応する 4 番目の係数です。
2 番目の主成分は縦軸にあり、変数 、 および に対する負の係数と、変数 に対する正の係数があります。
この 2 次元バイプロットには、13 件の各観測値に対応する点と、プロットの 2 つの主成分に対する各観測値のスコアを示す座標が含まれています。たとえば、プロットの左端に近い点は、第 1 主成分に対するスコアが最小になっています。各点は得点の最大値と係数の最大長を基準にスケーリングされるため、プロットからは各点の相対位置のみを決定できます。
ホテリングの T 二乗統計量を求めます。
標本データ セットを読み込みます。
load hald成分データには、4 つの変数に対しての 13 の観測値があります。
主成分分析を実行し、T 二乗値を要求します。
[coeff,score,latent,tsquared] = pca(ingredients); tsquared
tsquared = 13×1
5.6803
3.0758
6.0002
2.6198
3.3681
0.5668
3.4818
3.9794
2.6086
7.4818
4.1830
2.2327
2.7216
⋮
最初の 2 つの主成分のみを要求し、要求された主成分の縮小空間のみの T 二乗値を計算します。
[coeff,score,latent,tsquared] = pca(ingredients,'NumComponents',2);
tsquaredtsquared = 13×1
5.6803
3.0758
6.0002
2.6198
3.3681
0.5668
3.4818
3.9794
2.6086
7.4818
4.1830
2.2327
2.7216
⋮
主成分の縮小空間を指定する場合でも、pca は 4 つすべての主成分を使用して全体空間の T 二乗値を計算する点に注意してください。
縮小空間での T 二乗値は、その縮小空間におけるマハラノビス距離に対応します。
tsqreduced = mahal(score,score)
tsqreduced = 13×1
3.3179
2.0079
0.5874
1.7382
0.2955
0.4228
3.2457
2.6914
1.3619
2.9903
2.4371
1.3788
1.5251
⋮
全体空間の T 二乗値と縮小空間のマハラノビス距離の相違を計算することにより、破棄された空間において T 二乗値を計算します。
tsqdiscarded = tsquared - tsqreduced
tsqdiscarded = 13×1
2.3624
1.0679
5.4128
0.8816
3.0726
0.1440
0.2362
1.2880
1.2467
4.4915
1.7459
0.8539
1.1965
⋮
主成分によって説明される変化の比率を求めます。主成分空間におけるデータ表現を表示します。
標本データ セットを読み込みます。
load imports-85データ行列 X には、3 列目から 15 列目に 13 個の連続変数が含まれています (wheel-base、length、width、height、curb-weight、engine-size、bore、stroke、compression-ratio、horsepower、peak-rpm、city-mpg、highway-mpg)。
これら変数の主成分によって説明される変化の比率を求めます。
[coeff,score,latent,tsquared,explained] = pca(X(:,3:15)); explained
explained = 13×1
64.3429
35.4484
0.1550
0.0379
0.0078
0.0048
0.0013
0.0011
0.0005
0.0002
0.0002
0.0000
0.0000
⋮
最初の 3 つの成分によって変動性全体の 99.95% が説明されています。
最初の 3 つの主成分からなる空間におけるデータ表現を可視化します。
scatter3(score(:,1),score(:,2),score(:,3)) axis equal xlabel('1st Principal Component') ylabel('2nd Principal Component') zlabel('3rd Principal Component')
データから、最大の変動性は 1 番目の主成分軸に沿っていることがわかります。これは、1 番目の軸として可能なすべての選択肢の中で、可能な最大の分散です。2 番目の主成分軸に沿った変動性は、2 番目の軸として可能な残りのすべての選択肢の中で最大です。3 番目の主成分軸には、2 番目の主成分軸に沿った変動性より有意に小さい、3 番目に大きい変動性があります。4 番目から 13 番目の主成分軸は、データのすべての変動性の 0.05% しか説明しないので、調べる価値はありません。
出力をスキップするには、それに対応する要素で代わりに ~ を使用します。たとえば、T 二乗値を求めない場合は次のように指定します。
[coeff,score,latent,~,explained] = pca(X(:,3:15));
あるデータ セットの主成分を求め、別のデータ セットに PCA を適用します。この手順は、機械学習モデル用の学習データ セットとテスト データ セットがある場合に便利です。たとえば、PCA を使用して学習データ セットを前処理してから、モデルに学習をさせることができます。テスト データ セットを使用して学習済みモデルを検定するには、学習データから取得した PCA 変換をテスト データ セットに適用する必要があります。
この例では、C/C++ コードを生成する方法も示します。pca はコード生成をサポートするので、学習データ セットを使用して PCA を実行するコードを生成し、PCA をテスト データ セットに適用できます。その後、コードをデバイスに展開します。このワークフローでは学習データを渡さなければなりませんが、サイズが非常に大きい可能性があります。デバイスのメモリを節約するため、学習と予測を分離することができます。MATLAB® で pca を使用し、生成されたコードによりデバイスで PCA を新しいデータに適用します。
C/C++ コードの生成には MATLAB® Coder™ が必要です。
新しいデータへの PCA の適用
readtable を使用して、データ セットを table に読み込みます。データ セットは CreditRating_Historical.dat ファイル内にあり、過去の信用格付けデータが格納されています。
creditrating = readtable('CreditRating_Historical.dat');
creditrating(1:5,:)ans=5×8 table
ID WC_TA RE_TA EBIT_TA MVE_BVTD S_TA Industry Rating
_____ _____ _____ _______ ________ _____ ________ _______
62394 0.013 0.104 0.036 0.447 0.142 3 {'BB' }
48608 0.232 0.335 0.062 1.969 0.281 8 {'A' }
42444 0.311 0.367 0.074 1.935 0.366 1 {'A' }
48631 0.194 0.263 0.062 1.017 0.228 4 {'BBB'}
43768 0.121 0.413 0.057 3.647 0.466 12 {'AAA'}
最初の列は各観測値の ID、最後の列は格付けです。2 ～ 7 番目の列を予測子データとして、最後の列 (Rating) を応答として指定します。
X = table2array(creditrating(:,2:7)); Y = creditrating.Rating;
最初の 100 個の観測値をテスト データとして、残りを学習データとして使用します。
XTest = X(1:100,:); XTrain = X(101:end,:); YTest = Y(1:100); YTrain = Y(101:end);
学習データ セット XTrain の主成分を求めます。
[coeff,scoreTrain,~,~,explained,mu] = pca(XTrain);
このコードは、coeff、scoreTrain、explained および mu という 4 つの出力を返します。explained (説明される分散の合計の比率) を使用して、変動性の 95% 以上を説明するために必要な成分の個数を求めます。coeff (主成分係数) と mu (XTrain の推定平均) を使用して、PCA をテスト データ セットに適用します。モデルに学習をさせるときに、XTrain ではなく scoreTrain (主成分スコア) を使用します。
主成分によって説明される変動性の比率を表示します。
explained
explained = 6×1
58.2614
41.2606
0.3875
0.0632
0.0269
0.0005
最初の 2 つの成分によって変動性全体の 95% 以上が説明されています。変動性の 95% 以上を説明するために必要な成分の個数を求めます。
idx = find(cumsum(explained)>95,1)
idx = 2
最初の 2 つの成分を使用して、分類木に学習をさせます。
scoreTrain95 = scoreTrain(:,1:idx); mdl = fitctree(scoreTrain95,YTrain);
mdl は ClassificationTree モデルです。
テスト セットに対して学習済みモデルを使用するには、学習データ セットから取得した PCA を使用して、テスト データ セットを変換する必要があります。mu を XTest から減算し、coeff を乗算することにより、テスト データ セットの主成分スコアを取得します。最初の 2 つの成分のスコアのみが必要なので、最初の 2 つの係数 coeff(:,1:idx) を使用します。
scoreTest95 = (XTest-mu)*coeff(:,1:idx);
学習済みモデル mdl と変換済みテスト データ セット scoreTest を関数 predict に渡して、テスト セットの格付けを予測します。
YTest_predicted = predict(mdl,scoreTest95);
コードの生成
PCA をデータに適用し学習済みモデルを使用して格付けを予測するコードを生成します。C/C++ コードの生成には MATLAB® Coder™ が必要であることに注意してください。
saveLearnerForCoder を使用して、分類モデルを myMdl.mat というファイルに保存します。
saveLearnerForCoder(mdl,'myMdl');テスト データ セット (XTest) と PCA の情報 (coeff および mu) を受け入れ、テスト データの格付けを返す、myPCAPredict という名前のエントリポイント関数を定義します。
MATLAB のアルゴリズムについてのコードを生成しようとしていることを指示するため、コンパイラ命令 %#codegen (またはプラグマ) をエントリポイント関数のシグネチャの後に追加します。この命令を追加すると、コード生成時にエラーになる違反の診断と修正を MATLAB Code Analyzer が支援します。
function label = myPCAPredict(XTest,coeff,mu) %#codegen % Transform data using PCA scoreTest = bsxfun(@minus,XTest,mu)*coeff; % Load trained classification model mdl = loadLearnerForCoder('myMdl'); % Predict ratings using the loaded model label = predict(mdl,scoreTest);
myPCAPredict は、coeff と mu を使用して PCA を新しいデータに適用してから、変換済みデータを使用して格付けを予測します。この方法では、(サイズが非常に大きい可能性がある) 学習データを渡しません。
メモ: このページの右上にあるボタンをクリックしてこの例を MATLAB® で開くと、MATLAB® で例のフォルダーが開きます。このフォルダーには、エントリポイント関数のファイルが含まれています。
codegen (MATLAB Coder)を使用してコードを生成します。C および C++ は静的な型の言語なので、エントリポイント関数内のすべての変数のプロパティをコンパイル時に決定しなければなりません。データ型と正確な入力配列サイズを指定するため、-args オプションを使用して、特定のデータ型および配列サイズをもつ一連の値を表す MATLAB® 式を渡します。コンパイル時に観測値の個数が不明である場合、coder.typeof (MATLAB Coder) を使用して可変サイズの入力を指定することもできます。詳細については、コード生成用の可変サイズ引数の指定を参照してください。
codegen myPCAPredict -args {coder.typeof(XTest,[Inf,6],[1,0]),coeff(:,1:idx),mu}
Code generation successful.
codegen は、プラットフォームに依存する拡張子をもつ MEX 関数 myPCAPredict_mex を生成します。
生成されたコードを検証します。
YTest_predicted_mex = myPCAPredict_mex(XTest,coeff(:,1:idx),mu); isequal(YTest_predicted,YTest_predicted_mex)
ans = logical
1
isequal は、すべての入力が等しいことを意味する logical 1 (true) を返します。この比較により、mdl の関数 predict と関数 myPCAPredict_mex が同じ格付けを返すことを確認します。
コード生成の詳細については、コード生成の紹介および機械学習アプリからエクスポートされたモデルを使用したコマンド ラインでのコード生成を参照してください。後者では、分類学習器アプリを使用して PCA を実行しモデルに学習をさせる方法と、学習済みモデルに基づいて新しいデータのラベルを予測する C/C++ コードを生成する方法を説明しています。
入力引数
名前と値の引数
出力引数
詳細
アルゴリズム
関数 pca は、符号規則を適用して、coefs の各列で大きさが最大である要素を正にします。係数ベクトルの符号を変更しても意味は変わりません。
代替機能
アプリ
pca をライブ エディターで対話的に実行するには、[次元削減] ライブ エディター タスクを使用します。
参照
[1] Jolliffe, I. T. Principal Component Analysis. 2nd ed., Springer, 2002.
[2] Krzanowski, W. J. Principles of Multivariate Analysis. Oxford University Press, 1988.
[3] Seber, G. A. F. Multivariate Observations. Wiley, 1984.
[4] Jackson, J. E. A. User's Guide to Principal Components. Wiley, 1988.
[5] Roweis, S. “EM Algorithms for PCA and SPCA.” In Proceedings of the 1997 Conference on Advances in Neural Information Processing Systems. Vol.10 (NIPS 1997), Cambridge, MA, USA: MIT Press, 1998, pp. 626–632.
[6] Ilin, A., and T. Raiko. “Practical Approaches to Principal Component Analysis in the Presence of Missing Values.” J. Mach. Learn. Res.. Vol. 11, August 2010, pp. 1957–2000.
拡張機能
バージョン履歴
R2012b で導入
参考
barttest | biplot | canoncorr | factoran | pcacov | pcares | rotatefactors | ppca | incrementalPCA | 次元削減
