Ⅰ 鍝浣嶅ソ蹇冧漢鍙浠ユ彁渚涗竴浠絃DA浜鴻劯璇嗗埆鐨刴atlab紼嬪簭鍟婏紵璋㈣阿浜嗭紒
浠ヤ笅鏄疞DA鐨刴鏂囦歡鍑芥暟錛
浣犵◢紼嶆敼鏀瑰氨鑳界敤浜嗭紒
function [eigvector, eigvalue, elapse] = LDA(gnd,options,data)
% LDA: Linear Discriminant Analysis
%
% [eigvector, eigvalue] = LDA(gnd, options, data)
%
% Input:
% data - Data matrix. Each row vector of fea is a data point.
% gnd - Colunm vector of the label information for each
% data point.
% options - Struct value in Matlab. The fields in options
% that can be set:
%
% Regu - 1: regularized solution,
% a* = argmax (a'X'WXa)/(a'X'Xa+ReguAlpha*I)
% 0: solve the sinularity problem by SVD
% Default: 0
%
% ReguAlpha - The regularization parameter. Valid
% when Regu==1. Default value is 0.1.
%
% ReguType - 'Ridge': Tikhonov regularization
% 'Custom': User provided
% regularization matrix
% Default: 'Ridge'
% regularizerR - (nFea x nFea) regularization
% matrix which should be provided
% if ReguType is 'Custom'. nFea is
% the feature number of data
% matrix
% Fisherface - 1: Fisherface approach
% PCARatio = nSmp - nClass
% Default: 0
%
% PCARatio - The percentage of principal
% component kept in the PCA
% step. The percentage is
% calculated based on the
% eigenvalue. Default is 1
% (100%, all the non-zero
% eigenvalues will be kept.
% If PCARatio > 1, the PCA step
% will keep exactly PCARatio principle
% components (does not exceed the
% exact number of non-zero components).
%
%
% Output:
% eigvector - Each column is an embedding function, for a new
% data point (row vector) x, y = x*eigvector
% will be the embedding result of x.
% eigvalue - The sorted eigvalue of LDA eigen-problem.
% elapse - Time spent on different steps
%
% Examples:
%
% fea = rand(50,70);
% gnd = [ones(10,1);ones(15,1)*2;ones(10,1)*3;ones(15,1)*4];
% options = [];
% options.Fisherface = 1;
% [eigvector, eigvalue] = LDA(gnd, options, fea);
% Y = fea*eigvector;
%
%
% See also LPP, constructW, LGE
%
%
%
%Reference:
%
% P. N. Belhumeur, J. P. Hespanha, and D. J. Kriegman, 鎻坕genfaces
% vs. fisherfaces: recognition using class specific linear
% projection,� IEEE Transactions on Pattern Analysis and Machine
% Intelligence, vol. 19, no. 7, pp. 711-720, July 1997.
%
% Deng Cai, Xiaofei He, Yuxiao Hu, Jiawei Han, and Thomas Huang,
% "Learning a Spatially Smooth Subspace for Face Recognition", CVPR'2007
%
% Deng Cai, Xiaofei He, Jiawei Han, "SRDA: An Efficient Algorithm for
% Large Scale Discriminant Analysis", IEEE Transactions on Knowledge and
% Data Engineering, 2007.
%
% version 2.1 --June/2007
% version 2.0 --May/2007
% version 1.1 --Feb/2006
% version 1.0 --April/2004
%
% Written by Deng Cai (dengcai2 AT cs.uiuc.e)
%
if ~exist('data','var')
global data;
end
if (~exist('options','var'))
options = [];
end
if ~isfield(options,'Regu') | ~options.Regu
bPCA = 1;
if ~isfield(options,'PCARatio')
options.PCARatio = 1;
end
else
bPCA = 0;
if ~isfield(options,'ReguType')
options.ReguType = 'Ridge';
end
if ~isfield(options,'ReguAlpha')
options.ReguAlpha = 0.1;
end
end
tmp_T = cputime;
% ====== Initialization
[nSmp,nFea] = size(data);
if length(gnd) ~= nSmp
error('gnd and data mismatch!');
end
classLabel = unique(gnd);
nClass = length(classLabel);
Dim = nClass - 1;
if bPCA & isfield(options,'Fisherface') & options.Fisherface
options.PCARatio = nSmp - nClass;
end
if issparse(data)
data = full(data);
end
sampleMean = mean(data,1);
data = (data - repmat(sampleMean,nSmp,1));
bChol = 0;
if bPCA & (nSmp > nFea+1) & (options.PCARatio >= 1)
DPrime = data'*data;
DPrime = max(DPrime,DPrime');
[R,p] = chol(DPrime);
if p == 0
bPCA = 0;
bChol = 1;
end
end
%======================================
% SVD
%======================================
if bPCA
if nSmp > nFea
ddata = data'*data;
ddata = max(ddata,ddata');
[eigvector_PCA, eigvalue_PCA] = eig(ddata);
eigvalue_PCA = diag(eigvalue_PCA);
clear ddata;
maxEigValue = max(abs(eigvalue_PCA));
eigIdx = find(eigvalue_PCA/maxEigValue < 1e-12);
eigvalue_PCA(eigIdx) = [];
eigvector_PCA(:,eigIdx) = [];
[junk, index] = sort(-eigvalue_PCA);
eigvalue_PCA = eigvalue_PCA(index);
eigvector_PCA = eigvector_PCA(:, index);
%=======================================
if options.PCARatio > 1
idx = options.PCARatio;
if idx < length(eigvalue_PCA)
eigvalue_PCA = eigvalue_PCA(1:idx);
eigvector_PCA = eigvector_PCA(:,1:idx);
end
elseif options.PCARatio < 1
sumEig = sum(eigvalue_PCA);
sumEig = sumEig*options.PCARatio;
sumNow = 0;
for idx = 1:length(eigvalue_PCA)
sumNow = sumNow + eigvalue_PCA(idx);
if sumNow >= sumEig
break;
end
end
eigvalue_PCA = eigvalue_PCA(1:idx);
eigvector_PCA = eigvector_PCA(:,1:idx);
end
%=======================================
eigvalue_PCA = eigvalue_PCA.^-.5;
data = (data*eigvector_PCA).*repmat(eigvalue_PCA',nSmp,1);
else
ddata = data*data';
ddata = max(ddata,ddata');
[eigvector, eigvalue_PCA] = eig(ddata);
eigvalue_PCA = diag(eigvalue_PCA);
clear ddata;
maxEigValue = max(eigvalue_PCA);
eigIdx = find(eigvalue_PCA/maxEigValue < 1e-12);
eigvalue_PCA(eigIdx) = [];
eigvector(:,eigIdx) = [];
[junk, index] = sort(-eigvalue_PCA);
eigvalue_PCA = eigvalue_PCA(index);
eigvector = eigvector(:, index);
%=======================================
if options.PCARatio > 1
idx = options.PCARatio;
if idx < length(eigvalue_PCA)
eigvalue_PCA = eigvalue_PCA(1:idx);
eigvector = eigvector(:,1:idx);
end
elseif options.PCARatio < 1
sumEig = sum(eigvalue_PCA);
sumEig = sumEig*options.PCARatio;
sumNow = 0;
for idx = 1:length(eigvalue_PCA)
sumNow = sumNow + eigvalue_PCA(idx);
if sumNow >= sumEig
break;
end
end
eigvalue_PCA = eigvalue_PCA(1:idx);
eigvector = eigvector(:,1:idx);
end
%=======================================
eigvalue_PCA = eigvalue_PCA.^-.5;
eigvector_PCA = (data'*eigvector).*repmat(eigvalue_PCA',nFea,1);
data = eigvector;
clear eigvector;
end
else
if ~bChol
DPrime = data'*data;
% options.ReguAlpha = nSmp*options.ReguAlpha;
switch lower(options.ReguType)
case {lower('Ridge')}
for i=1:size(DPrime,1)
DPrime(i,i) = DPrime(i,i) + options.ReguAlpha;
end
case {lower('Tensor')}
DPrime = DPrime + options.ReguAlpha*options.regularizerR;
case {lower('Custom')}
DPrime = DPrime + options.ReguAlpha*options.regularizerR;
otherwise
error('ReguType does not exist!');
end
DPrime = max(DPrime,DPrime');
end
end
[nSmp,nFea] = size(data);
Hb = zeros(nClass,nFea);
for i = 1:nClass,
index = find(gnd==classLabel(i));
classMean = mean(data(index,:),1);
Hb (i,:) = sqrt(length(index))*classMean;
end
elapse.timeW = 0;
elapse.timePCA = cputime - tmp_T;
tmp_T = cputime;
if bPCA
[mpVec,eigvalue,eigvector] = svd(Hb,'econ');
eigvalue = diag(eigvalue);
eigIdx = find(eigvalue < 1e-3);
eigvalue(eigIdx) = [];
eigvector(:,eigIdx) = [];
eigvalue = eigvalue.^2;
eigvector = eigvector_PCA*(repmat(eigvalue_PCA,1,length(eigvalue)).*eigvector);
else
WPrime = Hb'*Hb;
WPrime = max(WPrime,WPrime');
dimMatrix = size(WPrime,2);
if Dim > dimMatrix
Dim = dimMatrix;
end
if isfield(options,'bEigs')
if options.bEigs
bEigs = 1;
else
bEigs = 0;
end
else
if (dimMatrix > 1000 & Dim < dimMatrix/10) | (dimMatrix > 500 & Dim < dimMatrix/20) | (dimMatrix > 250 & Dim < dimMatrix/30)
bEigs = 1;
else
bEigs = 0;
end
end
if bEigs
%disp('use eigs to speed up!');
option = struct('disp',0);
if bChol
option.cholB = 1;
[eigvector, eigvalue] = eigs(WPrime,R,Dim,'la',option);
else
[eigvector, eigvalue] = eigs(WPrime,DPrime,Dim,'la',option);
end
eigvalue = diag(eigvalue);
else
[eigvector, eigvalue] = eig(WPrime,DPrime);
eigvalue = diag(eigvalue);
[junk, index] = sort(-eigvalue);
eigvalue = eigvalue(index);
eigvector = eigvector(:,index);
if Dim < size(eigvector,2)
eigvector = eigvector(:, 1:Dim);
eigvalue = eigvalue(1:Dim);
end
end
end
for i = 1:size(eigvector,2)
eigvector(:,i) = eigvector(:,i)./norm(eigvector(:,i));
end
elapse.timeMethod = cputime - tmp_T;
elapse.timeAll = elapse.timePCA + elapse.timeMethod;
Ⅱ matlab中的降維函數是什麼
drttoolbox : Matlab Toolbox for Dimensionality Rection是Laurens van der Maaten數據降維的工具箱。
裡面囊括了幾乎所有的數據降維演算法：
- Principal Component Analysis ('PCA')
- Linear Discriminant Analysis ('LDA')
- Independent Component Analysis ('ICA')
- Multidimensional scaling ('MDS')
- Isomap ('Isomap')
- Landmark Isomap ('LandmarkIsomap')
- Locally Linear Embedding ('LLE')
- Locally Linear Coordination ('LLC')
- Laplacian Eigenmaps ('Laplacian')
- Hessian LLE ('HessianLLE')
- Local Tangent Space Alignment ('LTSA')
- Diffusion maps ('DiffusionMaps')
- Kernel PCA ('KernelPCA')
- Generalized Discriminant Analysis ('KernelLDA')
- Stochastic Neighbor Embedding ('SNE')
- Neighborhood Preserving Embedding ('NPE')
- Linearity Preserving Projection ('LPP')
- Stochastic Proximity Embedding ('SPE')
- Linear Local Tangent Space Alignment ('LLTSA')
- Simple PCA ('SPCA')
Ⅲ 常用降維方法之PCA 和 LDA
PCA本質上是將方差最大的方向作為主要特徵，並且在各個正交方向上將數據「離相關」，也就是讓它們在不同正交方向上沒有相關性。而方差最大的那個維度是主成分。
PCA是比較常見的線性降維方法,通過線性投影將高維數據映射到低維數據中,所期望的是在投影的維度上,新特徵自身的方差盡量大,方差越大特徵越有效,盡量使產生的新特徵間的相關性越小。
PCA演算法的具體操作為對所有的樣本進行中心化操作,計算樣本的協方差矩陣,然後對協方差矩陣做特徵值分解,取最大的n個特徵值對應的特徵向量構造投影矩陣。
再舉個栗子：
下面舉一個簡單的例子，說明PCA的過程。
假設我們的數據集有10個二維數據(2.5,2.4), (0.5,0.7), (2.2,2.9), (1.9,2.2), (3.1,3.0), (2.3, 2.7), (2, 1.6), (1, 1.1), (1.5, 1.6), (1.1, 0.9)，需要用PCA降到1維特徵。
首先我們對樣本中心化，這里樣本的均值為(1.81, 1.91),所有的樣本減去這個均值向量後，即中心化後的數據集為(0.69, 0.49), (-1.31, -1.21), (0.39, 0.99), (0.09, 0.29), (1.29, 1.09), (0.49, 0.79), (0.19, -0.31), (-0.81, -0.81), (-0.31, -0.31), (-0.71, -1.01)。
現在我們開始求樣本的協方差矩陣，由於我們是二維的，則協方差矩陣為：
對於我們的數據，求出協方差矩陣為：
求出特徵值為（0.0490833989， 1.28402771），對應的特徵向量分別為：
由於最大的k=1個特徵值為1.28402771，對於的k=1個特徵向量為 則我們的W=
我們對所有的數據集進行投影 得到PCA降維後的10個一維數據集為：(-0.827970186， 1.77758033， -0.992197494， -0.274210416， -1.67580142， -0.912949103， 0.0991094375， 1.14457216, 0.438046137， 1.22382056)
在上面的PCA演算法中，我們假設存在一個線性的超平面，可以讓我們對數據進行投影。但是有些時候，數據不是線性的，不能直接進行PCA降維。這里就需要用到和支持向量機一樣的核函數的思想，先把數據集從n維映射到線性可分的高維N>n,然後再從N維降維到一個低維度n', 這里的維度之間滿足n'<n<N。
使用了核函數的主成分分析一般稱之為核主成分分析(Kernelized PCA, 以下簡稱KPCA。假設高維空間的數據是由n維空間的數據通過映射ϕ產生。
則對於n維空間的特徵分解：
映射為：
通過在高維空間進行協方差矩陣的特徵值分解，然後用和PCA一樣的方法進行降維。一般來說，映射ϕ不用顯式的計算，而是在需要計算的時候通過核函數完成。由於KPCA需要核函數的運算，因此它的計算量要比PCA大很多。
這里對PCA演算法做一個總結。作為一個非監督學習的降維方法，它只需要特徵值分解，就可以對數據進行壓縮，去噪。因此在實際場景應用很廣泛。為了克服PCA的一些缺點，出現了很多PCA的變種，比如第六節的為解決非線性降維的KPCA，還有解決內存限制的增量PCA方法Incremental PCA，以及解決稀疏數據降維的PCA方法Sparse PCA等。
PCA演算法的主要優點有：
LDA（線性判別分析，Linear Discriminant Analysis）是另一種常用的降維方法，它是有監督的。LDA在模式識別領域（比如人臉識別，艦艇識別等圖形圖像識別領域）中有非常廣泛的應用，因此我們有必要了解下它的演算法原理。這里需要注意的是，此處的LDA與文本主題模型中的LDA（隱含狄利克雷分布，Latent Dirichlet Allocation）並不相同，他是一種處理文檔的主題模型。
LDA是一種監督學習的降維技術，也就是說它的數據集的每個樣本是有類別輸出的。這點和PCA不同。PCA是不考慮樣本類別輸出的無監督降維技術。
LDA的思想可以用一句話概括，就是「投影後類內方差最小，類間方差最大」。
什麼意思呢？ 我們要將數據在低維度上進行投影，投影後希望每一種類別數據的投影點盡可能的接近，而不同類別的數據的類別中心之間的距離盡可能的大。
可能還是有點抽象，我們先看看最簡單的情況。假設我們有兩類數據 分別為紅色和藍色，如下圖所示，這些數據特徵是二維的，我們希望將這些數據投影到一維的一條直線，讓每一種類別數據的投影點盡可能的接近，而紅色和藍色數據中心之間的距離盡可能的大。
以上就是使用LDA進行降維的演算法流程。實際上LDA除了可以用於降維以外，還可以用於分類。一個常見的LDA分類基本思想是假設各個類別的樣本數據符合高斯分布，這樣利用LDA進行投影後，可以利用極大似然估計計算各個類別投影數據的均值和方差，進而得到該類別高斯分布的概率密度函數。當一個新的樣本到來後，我們可以將它投影，然後將投影後的樣本特徵分別帶入各個類別的高斯分布概率密度函數，計算它屬於這個類別的概率，最大的概率對應的類別即為預測類別。
LDA用於降維，和PCA有很多相同，也有很多不同的地方，因此值得好好的比較一下兩者的降維異同點。
這點可以從下圖形象的看出，在某些數據分布下LDA比PCA降維較優。
當然，某些某些數據分布下PCA比LDA降維較優，如下圖所示：
LDA演算法既可以用來降維，又可以用來分類，但是目前來說，主要還是用於降維。在我們進行圖像識別圖像識別相關的數據分析時，LDA是一個有力的工具。下面總結下LDA演算法的優缺點。
LDA演算法的主要優點有：
參考文章： 劉建平老師的博客園
Ⅳ 降維演算法之LDA(線性判別降維演算法)--有監督
LDA在模式識別領域( 比如人臉識別，艦艇識別等圖形圖像識別領域 )中有非常廣泛的應用，因此我們有必要了解下它的演算法原理。
不同於PCA方差最大化理論， LDA演算法的思想是將數據投影到低維空間之後，使得同一類數據盡可能的緊湊，不同類的數據盡可能的分散 。因此，LDA演算法是一種有監督的機器學習演算法。同時，LDA有如下兩個假設：(1)原始數據根據樣本均值進行分類。(2)不同類的數據擁有相同的協方差矩陣。當然，在實際情況中，不可能滿足以上兩個假設。但是 當數據主要是由均值來區分的時候，LDA一般都可以取得很好的效果 。
（1）計算類內散度矩陣
（2）計算類間散度矩陣
（3）計算矩陣
（4）對矩陣 進行特徵分解，計算最大的d個最大的特徵值對應的特徵向量組成W。
（5）計算投影後的數據點
以上就是使用LDA進行降維的演算法流程。實際上LDA除了可以用於降維以外，還可以用於分類。 一個常見的LDA分類基本思想是假設各個類別的樣本數據符合高斯分布 ， 這樣利用LDA進行投影後，可以利用極大似然估計計算各個累唄投影數據的均值和方差，進而得到該類別高斯分布的概率密度函數 。當一個新的樣本到來後，我們可以將它投影，然後將投影後的樣本特徵分別帶入各個類別的高斯分布概率密度函數，計算它屬於這個類別的概率，最大的概率對應的類別即為預測類別。LDA應用於分類現在似乎也不是那麼流行。
class sklearn.discriminant_analysis.LinearDiscriminantAnalysis(solver='svd', shrinkage=None, priors=None, n_components=None, store_covariance=False, tol=0.0001)
參數：
（1）solver: str類型，默認值為"svd"，
svd:使用奇異值分解求解，不用計算協方差矩陣，適用於特徵數量很大的情形，無法使用參數收縮(shrinkage)。
lsqr:最小平方QR分解，可以結合shrinkage使用。
eigen：特徵值分解，可以結合shrinkage使用。
（2）shrinkage: str or float類型，默認值為None
是否使用參數收縮
None:不使用參數收縮
auto:str，使用Ledoit-Wolf lemma
浮點數：自定義收縮比例。
（3）components：int類型，需要保留的特徵個數，小於等於n-1
屬性：
（1）covariances_：每個類的協方差矩陣，shape = [n_features, n_features]
（2）means_：類均值，shape = [n_features, n_feateures]
（3）priors_：歸一化的先驗概率。
（4）rotations_：LDA分析得到的主軸，shape = [n_features, n_component]
（5）scalings_：數組列表，每個高斯分布的方差σ
特點：
降維之後的維數最多為類別數-1。所以當數據維度很高，但是類別數少的時候，演算法並不適用 。LDA演算法既可以用來降維，又可以用來分類。但是目前來說，主要還是用於降維。在我們 進行圖像識別相關的數據分析時，LDA是一個有力的工具 。
優點：
（1） LDA在樣本分類信息依賴均值而不是方差的時候，比PCA之類的演算法較優 。
（2）在降維過程中可以使用類別的先驗知識經驗，而像PCA這樣的無監督學習則無法使用類別先驗知識。
缺點：
（1）LDA不適合非高斯分布樣本進行降維，PCA也存在這個問題。
（2）LDA降維最多降到類別數K-1的維數，如果我們降維的維度大於k-1，則不能使用LDA。 當然目前有一些LDA的進化版演算法可以繞過這個問題 。
（3） LDA在樣本分類信息依賴方差而不是均值的時候，降維效果不好 。
（4）LDA可能過度擬合數據。
二者都有 降維 的作用。
1.左 邊是PCA，屬於無監督方法 ，當數據沒有標簽時可以用它。 右邊是LDA，屬於監督學習方法 。考慮了數據的分類信息，這樣數據在低維空間上就可以分類了，減少了很多的運算量。
2. PCA主要是從特徵的協方差角度考慮，追求的是在降維之後能夠最大化保持數據的內在信息 。它不考慮分類信息，因此降低維度後，信息損失降到最低，但分類上可能會變得更加困難。 LDA追求的是降維後的數據點盡可能容易被區分 。降維後的樣本數據在新的維度空間有最大的類間距離和最小的類內方差，數據在低維空間有最佳的可分離性。
3. PCA降維後的維度數目是和數據維度相關的 ，原始數據是n維，那麼PCA後維度為1、2~n維。 LDA後的維度數目是和類別的個數相關的 ，原始數據是n維，一共有C個類別，那麼LDA後維度為1、2~C-1維。
4. PCA投影的坐標系都是正交的 。 LDA關注分類能力，不保證投影到的坐標系是正交的 。