PCA介绍
主成分分析 (Principal Component Analysis ,简称 PCA)是最常用的一种降维方法。
算法介绍
代码
function [p_d_mat,project_mat,allmean]=PCA2(initialize_mat,p)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%输入:
% initialize_mat: 原高维矩阵 m*n m代表特征 n代表样本数
% p: 要求降维后矩阵的维数
%输出:
%p_d_mat: 降维后的p维矩阵
%project_mat: 投影矩阵
%allmean: 均值
%%%%%%%%%%%%%%%%%%%%%%%% by LiSR %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%关于PCA2(initialize_mat,p)处理函数的说明
%1、对initialize_mat每一个特征维度求均值,然后减去均值得meanmat
%2、对meanmat'*meanmat求特征向量与特征值
%3、找出最大的特征值与特征向量,计算投影矩阵,project_mat=meanmat*max_v/sqrt(sigma) (为了特征向量单位阵)
%4、计算降维后的矩阵 p_d_mat=project_mat'*meanmat;
%5、除以p_d_mat矩阵最大的奇异值
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[features,rows]=size(initialize_mat);%原矩阵的特征(维)数,和列数 size返回行与列
allmean=mean(initialize_mat,2);%按列求原矩阵的均值 2代表按列求均值
meanmat=zeros(features,rows);%初始化均值矩阵
%计算差值,原矩阵所有列减去均值
for i=1:rows
meanmat(:,i)=initialize_mat(:,i)-allmean;
end
%利用差值矩阵构建协方差矩阵,即C=meanmat*meanmat';
%然后计算C的特征值和特征向量。由于C的维数过大,计算代价太大,
%利用奇异值分解原理求得协方差矩阵C=meanmat*meanmat'的特征值和特征向量
%即,利用求meanmat'*meanmat的特征值、向量求原协方差矩阵meanmat*meanmat'
%的特征值和特征向量
%[U,S,V]=svd(meanmat)
%注:Matlab中,求解出来的特征值是默认从大到小排序了的,而自然的,特征向量也是与排序后特征值相对应的。
[v,smat]=eig(meanmat'*meanmat);%求得特征向量所组成的v矩阵和特征值
%把特征值取出,放到数组中,按升序排列
w=1:rows;
for i=1:rows
w(i)=smat(i,i);
end
project_mat=zeros(features,p);%初始化投影矩阵
max_v=v(:,rows-p+1:rows); %取最大的特征值对应的特征向量
%计算奇异值
sigma=sqrt(w(rows-p+1:rows));
%计算投影矩阵,U_i=A*V_i/sqrt(sigma) 对其归一化
for i=1:p
project_mat(:,i)=meanmat*max_v(:,i)/sigma(i);
end
%计算降维后的矩阵
p_d_mat=project_mat'*meanmat;
%L_2 归一 norm->返回A的最大奇异值,即max(svd(A))
for i=1:size(p_d_mat,2)
p_d_mat(:,i)=p_d_mat(:,i)/norm(p_d_mat(:,i)); %假设A是一个矩阵,那么norm(A)或者norm(A,2)计算的就是A的2范数;
end
end
主函数
%pca+svm测试
%% a litte clean work
close all;
clear;
clc;
format compact;
%% load data 加载PaviaU数据集
%数据说明:训练集共有1800个样本,测试集有900个样本
%数据样本共有9种地表类型
load('PaviaU_gt.mat');
load('PaviaU.mat');
clsCnt = 9;
trnNum = 200;
%转化成svm训练所需要的格式
train_x = [];
train_y = [];
test_x = [];
test_y = [];
mm = max(max(max(paviaU)));
for i = 1 :clsCnt
temp = getpixel(paviaU,paviaU_gt,i);
clsNum = size(temp,1);
%trnNum = ceil(clsNum*trnPer);
tstNum = clsNum - trnNum;
train_x = [ train_x ; temp(1:trnNum,:) ];
cls_lab = zeros(trnNum,1);
cls_lab(:,1) = i;
train_y = [ train_y ; cls_lab];
test_x = [ test_x ; temp(trnNum+1:trnNum+100,:) ];
cls_lab = zeros(100,1);
cls_lab(:,1) = i;
test_y = [ test_y ; cls_lab];
end
%首先求取特征数
feature_num=size(train_x,2);
%准确度记录
Accuracy_record_train = [];
Accuracy_record_test = [];
%从103个特征中选取1~103个特征
for num=1:20
tic
%分别是训练集与测试集以及相应的标签
traindata = double(train_x/mm);
testdata = double(test_x/mm);
trainlabel=train_y;
testlabel=test_y;
%% 进行pca数据预处理
fprintf('iter %d',num);
%两种计算PCA的方法
[traindata_pca,project_mat,allmean]=PCA_mine(traindata',num);
test_samples=size(testdata,1);
for i=1:test_samples
testdata_meanmat(:,i)=testdata(i,:)-allmean'; %testdata 103*1
end
testdata_pca=project_mat'*testdata_meanmat;
% 利用训练集合建立分类模型
%-s 代表svm类型
%-t代表 核函数类型
%-c 设置C-SVC epsilon-SVR n-SVR惩罚参数 默认为1
%-g 设置核函数中 gamma的值
%model = svmtrain(trainlabel,traindata,'-s 0 -t 2 -c 1.2 -g 2.8');
%model = svmtrain(trainlabel,traindata,'-s 0 -t 1 -c 1.2 -g 2.8');
model = svmtrain(trainlabel,traindata_pca','-s 0 -t 1 -c 1.2 -g 2.8');
% 利用建立的svm模型看其在训练集合上的分类效果
[ptrain,acctrain,dec_values1] = svmpredict(trainlabel,traindata_pca',model);
% 预测测试集合标签
[ptest,acctest,dec_values2] = svmpredict(testlabel,testdata_pca',model);
%记录准确度
Accuracy_record_train=[Accuracy_record_train,acctrain(1)];
Accuracy_record_test=[Accuracy_record_test,acctest(1)];
%%
toc;
end
%% 绘制测试集 训练集准确度曲线
x = 1:1:20;
%x = 1:1:feature_num;
%y = sin(x);
plot(x,Accuracy_record_train)
hold on
y2 = cos(x);
plot(x,Accuracy_record_test,'r:')
grid on
xlabel('pca降维后的特征数')
ylabel('准确度')
%axis([xmin xmax ymin ymax]); % 设置坐标轴在指定的区间
%legend('训练集准确度','测试集准确度') %为图片添加图例
title('pca降维后的特征数与svm分类准确度的关系示意图') %添加图像标题
text(1.5,0.3,'测试集准确度') %将'测试集准确度'这个注解加到坐标中的某个位置
gtext('测试集准确度') % 用鼠标的光标定位,将'测试集准确度'这个注解放在你鼠标点击的地方
text(1.5,0.3,'训练集准确度') %将'测试集准确度'这个注解加到坐标中的某个位置
gtext('训练集准确度') % 用鼠标的光标定位,将'测试集准确度'这个注解放在你鼠标点击的地方