【深圳大學機器學習】實驗一：PCA算法

實驗目的

1、實現PCA算法的人臉重構，即用20,40,60,80,...,160個投影來重構圖像的效果。

2、實現PCA算法的人臉識別，給出不少于三個人臉數據庫上 10,20,30,...,160維的人臉識別識別率，并打印出識別率變化曲線圖。

3、用PCA用來進行人臉圖像降維，實現3個不同數據集多個子集（或多個類）的二維和三維空間實現數據的可視化，不同類用不同顏色表示。

4、比較PCA、2DPCA、B2PCA在各個人臉數據集上的最高識別率，比較重構圖像的差異。

實驗步驟：

一，基于PCA算法的人臉重構與人臉識別，以及數據可視化

1，（條件：基于AR數據集的前120人，每人取前三張圖像作為訓練數據，后七張作為識別數據）

原圖

不同維度下的重構效果（AR001-1.tif）

維度	正確率
20	62.70%
40	68.52%
60	71.33%
80	74.27%
100	74.53%
120	75.20%
140	75.41%
160	75.89%

取用不同數量特征值和正確率??????? 基于此情況下人臉識別識別率曲線

數據可視化

???????????? 二維空間?????????????????????? 三維空間

2，（條件：基于orl數據集的前40人，每人取前三張圖像作為訓練數據，后七張作為識別數據）

原圖

不同維度下重構效果（orl1_1.bmp）

維度	正確率
20	78.92%
40	83.57%
60	85.00%
80	85.71%
100	86.42%
120	86.07%
140	86.07%
160	86.07%

取用不同數量特征值和正確率??????? 基于此情況下人臉識別識別率曲線?????? ???

數據可視化

二維空間?????????????????????? 三維空間

3，（條件：基于Yale數據集的前15人，每人取前三張圖像作為訓練數據，后七張作為識別數據）

原圖

不同維度下重構效果(subject01_1.bmp)

維度	正確率
20	87.73%
40	86.77%
60	86.77%
80	86.77%
100	86.77%
120	86.77%
140	86.77%
160	86.77%

取用不同數量特征值和正確率???????? 基于此情況下人臉識別識別率曲線

數據可視化

???????????? 二維空間?????????????????????? 三維空間

二，PCA,2DPCA,B2DPCA重構圖像與識別率差異

條件：在ORL數據集上的重構效果（40個人，取前六張訓練，后四張識別的情況下展示重構效果）

重構效果：

原圖

PCA：

2DPCA:

B2DPCA

識別正確率（orl數據集，40個人，前六張訓練，后四張識別）：

PCA:

最高正確率為96.3%

2DPCA:

最高正確率為97.5%

B2DPCA:

最高正確率為97.5%

嘗試不同訓練圖像與識別圖像比例后結果如下：

訓練情況	訓練3，識別7	訓練4，識別6	訓練5，識別5	訓練6，識別4	訓練7，識別3
PCA	86.4%	89.6%	90.0%	96.3%	96.7%
2DPCA	93.2%	93.3%	93.5%	97.5%	97.5%
B2DPCA	93.2%	93.3%	93.5%	97.5%	97.5%

從表中可以知道，2DPCA和B2DPCA在ORL訓練集上的最高正確率基本相同，高于PCA的正確率，同時隨著訓練數據的增加，PCA與它們識別正確率的差距逐漸減小。

進一步發現，在訓練圖像6張，識別圖像四張的情況下，PCA，2DPCA,B2DPCA識別時間如下：

PCA	2DPCA	B2DPCA
4.09s	1.59s	1.22s

可以看出，2DPCA和B2DPCA雖然在正確率上差別不大，但B2DPCA卻比2DPCA更加節省時間，二者所耗時間均遠遠小于傳統PCA算法。

代碼：

PCA:

%1，測試與訓練圖像的錄入（附加計時以查看效率）
tic
reshaped_faces=[];
%訓練圖像讀取（以orl為例）
for i=1:40for j=1:10 a=imread(strcat('C:\Users\18023\Desktop\粵海小清華\homework\計導文件\人臉數據庫\ORL56_46\orl',num2str(i),'_',num2str(j),'.bmp'));     b = reshape(a,2576,1); %將每一張人臉拉成列向量b=double(b); %uint-8轉換為double類型，避免人臉展示時全灰的影響       reshaped_faces=[reshaped_faces, b];  %儲存每一列轉化后的列向量到矩陣reshaped_faces中end
end%用每個人前六張作為訓練圖，后四張作為測試圖
%初始化測試集集矩陣
test_data_index = [];
%初始化訓練集矩陣
train_data_index = [];  
for k=0:(size(reshaped_faces,2)/j)-1train_data_index= [train_data_index j*k+1:j*k+3];%建立一個用于查找所需訓練圖像的目錄test_data_index= [test_data_index j*k+4:j*(k+1)];%同理，是測試圖像的目錄
end
train_data = reshaped_faces(:,train_data_index);%根據index讀取矩陣 reshaped_faces的數據添加到矩陣train_data中
test_data = reshaped_faces(:, test_data_index);%根據index1讀取矩陣 reshaped_faces的數據添加到矩陣test_data中% 2-特征向量矩陣（all_eigen_face)求得過程
% 求平均臉
mean_face = mean(train_data,2);%用mean函數求訓練集矩陣的平均值,2是對每一行求平均，最后得到一個和原先尺寸一樣的一個人臉，每個像素是平均顏色
% 中心化
normalizeddata= (train_data - mean_face);%得到中心化人臉（此處減去一列等于traindata每一列都減去這一列） 
%PCA特征值分解的矩陣為去中心化的矩陣乘其轉置矩陣
cov= normalizeddata * normalizeddata';
%特征值分解
[eigen_vectors, eigenvalues] = eig(cov);
% 從對角矩陣獲取特征值
eigen_values = diag(eigenvalues);
% 對特征值按索引進行從大到小排序
[sorted_eigen_values, index] = sort(eigen_values, 'descend'); %此處建立index獲取特征值的位置，以便尋找對應的特征向量
%按照特征值的大小排列特征向量
sorted_eigen_vectors = eigen_vectors(:, index);
all_eigen_faces = sorted_eigen_vectors;%2重構過程
single_face = normalizeddata(:,1);%代碼中取第一張臉進行重構
index = 1;
figure;
for dim=[20,40,60,80,100,120,140,160]eigen_faces = all_eigen_faces(:,1:dim);%降維與再次重構rebuild_face = eigen_faces * (eigen_faces' * single_face) + mean_face;subplot(2,4, index); index = index + 1;fig = imshow(mat2gray(reshape(rebuild_face, [56, 46])));title(sprintf("維度=%d", dim));    if (dim == 160)waitfor(fig);end
end% 3.人臉識別
index = 1;       
k = length(10:10:160);  %選取適當間隔
Y = zeros(1, k); 
% 遍歷每個維度
for i = 1:kd = 10 * i;  % 通過i確定當前維度% 取出相應數量特征向量eigen_faces = all_eigen_faces(:, 1:d);% 測試、訓練數據降維projected_train_data = eigen_faces' * (train_data - mean_face);projected_test_data = eigen_faces' * (test_data - mean_face);% 測試前準備過程correct_predict_number = 0;test_face_number = size(projected_test_data, 2);% 讀取測試數據與初始化for each_test_face_index = 1:test_face_numbereach_test_face = projected_test_data(:, each_test_face_index);% 初始化最小距離min_distance = inf;predicted_label = -1;%在歐氏距離下與每一個降維后的數據進行比較for each_train_face_index = 1:size(projected_train_data, 2)% 計算距離distance = norm(each_test_face - projected_train_data(:, each_train_face_index));% 如果當前距離小于之前的最小距離，則更新最小距離和標簽（1NN，取最近的一個）if distance < min_distancemin_distance = distance;predicted_label = floor((train_data_index(1, each_train_face_index) - 1) / j) + 1;endend% 獲取該數據真實的標簽real_label = floor((test_data_index(1, each_test_face_index) - 1) / j) + 1;% 判斷預測結果是否正確if predicted_label == real_labelcorrect_predict_number = correct_predict_number + 1; % 如果預測正確，增加正確數endend % 計算正確率correct_rate = correct_predict_number / test_face_number;% 將當前維度對應的正確率存入YY(i) = correct_rate;  % 使用索引i來存儲正確率，確保Y與d一致% 輸出當前實驗結果fprintf("維度=%d，總測試樣本：%d，正確數:%d，正確率：%1f\n", d, test_face_number, correct_predict_number, correct_rate);
end
toc% 4，制作折線圖使正確率隨維度變化可視化
figure;
plot(10:10:160, Y, '-o');  % x軸為維度數，y軸為正確率
xlabel('維度數');
ylabel('正確率');
title('人臉識別正確率與維度變化');
waitfor(fig);%5，散點圖繪制
for t=[2 3]eigen_faces=all_eigen_faces(:,1:t);projected_test_data=eigen_faces'*(test_data - mean_face);color=[];for x=1:size(test_data,2)color=[color floor((x-1)/19*2)];endif (t == 2)waitfor(scatter(projected_test_data(1, :), projected_test_data(2, :), [], color));elsewaitfor(scatter3(projected_test_data(1, :), projected_test_data(2, :), projected_test_data(3, :), [], color));endend

2DPCA

%1，測試與訓練圖像的錄入
tic
reshaped_faces=[];
%訓練圖像讀取（以orl為例）
for i=1:40for j=1:10%a=imread(strcat('C:\Users\18023\Desktop\粵海小清華\homework\計導文件\人臉數據庫\ORL56_46\orl',num2str(i),'_',num2str(j),'.bmp'));     b = reshape(a,2576,1); b=double(b);       reshaped_faces=[reshaped_faces, b]; end
end%用每個人前六張作為訓練圖，后四張作為測試圖
%定義測試集集矩陣
test_data_index = [];
%定義訓練集矩陣
train_data_index = [];  
for k=0:(size(reshaped_faces,2)/j)-1train_data_index= [train_data_index j*k+1:j*k+6];%建立一個用于查找所需訓練圖像的目錄test_data_index= [test_data_index j*k+7:j*(k+1)];%同理，是測試圖像的目錄
end
train_data = reshaped_faces(:,train_data_index);%讀取矩陣 reshaped_faces的數據添加到矩陣train_data中
test_data = reshaped_faces(:, test_data_index);%讀取矩陣 reshaped_faces的數據添加到矩陣test_data中
mean_face = mean(train_data,2);%2，2DPCA矩陣求取
%對于訓練矩陣重新還原成二維疊加的狀態
train_matrices = cell(size(train_data, 2), 1);
for i = 1:size(train_data, 2)train_matrices{i} = reshape(train_data(:, i), [56, 46]);
end
%所有訓練圖像矩陣的均值
mean_face_matrix = zeros(56, 46);
for i = 1:length(train_matrices)mean_face_matrix = mean_face_matrix + train_matrices{i};
end
mean_face_matrix = mean_face_matrix / length(train_matrices);% 初始化2DPCA協方差矩陣
cov_2dpca = zeros(46, 46);  
% 計算2DPCA的協方差矩陣
for i = 1:length(train_matrices)centered_matrix = train_matrices{i} - mean_face_matrix;%每一個二維數據進行去中心化cov_contribution = centered_matrix' * centered_matrix;%求取其協方差矩陣cov_2dpca = cov_2dpca + cov_contribution;%每一個二維數據的協方差矩陣，都加到最后求解的矩陣上
end
cov_2dpca = cov_2dpca / length(train_matrices);%除以所有參與的二維矩陣數，求得均值
%cov_2dpca即為最終進行協方差分解的矩陣
%特征值分解與按照特征值從大到小排列
[eigen_vectors, eigenvalues] = eig(cov_2dpca);
eigen_values = diag(eigenvalues);
[sorted_eigen_values, index] = sort(eigen_values, 'descend');
sorted_eigen_vectors = eigen_vectors(:, index);
all_eigen_faces = sorted_eigen_vectors;figure
% 3，計算重構過程，展示第一張訓練圖像的重構
for idx = 1:4dim = idx * 10;eigen_faces = all_eigen_faces(:, 1:dim);  % 選擇前dim個特征向量% 計算投影系數test_image = train_matrices{1};  % 獲取第一張訓練圖像centered_image = test_image - mean_face_matrix;  % 中心化圖像projection = centered_image * eigen_faces;  % 投影到特征空間% 重構圖像rebuild_face = mean_face_matrix + (projection * eigen_faces');  % 重構% 顯示重構結果subplot(1, 4, idx);imshow(mat2gray(rebuild_face));  % 顯示重構圖像title(sprintf('維度=%d', dim));
end% 4.人臉識別
index = 1;       
k = length(10:10:40);
Y = zeros(1, k);
% 將測試圖像轉換為矩陣形式
test_matrices = cell(size(test_data, 2), 1);
for i = 1:size(test_data, 2)test_matrices{i} = reshape(test_data(:, i), [56, 46]);
end
% 遍歷每個維度
for i = 1:kd = 10* i;% 取出相應數量特征向量projection_vectors = all_eigen_faces(:, 1:d);% 計算訓練和測試圖像的特征train_features = zeros(size(train_data, 2), 56*d);test_features = zeros(size(test_data, 2), 56*d);% 對每張訓練圖像進行投影for j = 1:length(train_matrices)% 計算投影projected = train_matrices{j} * projection_vectors; % 56×d% 將投影結果展平為特征向量train_features(j, :) = projected(:)';end% 對每張測試圖像進行投影for j = 1:length(test_matrices)projected = test_matrices{j} * projection_vectors;test_features(j, :) = projected(:)';end% 測試前準備過程correct_predict_number = 0;test_face_number = size(test_features, 1);% 遍歷每一個待測試人臉for each_test_face_index = 1:test_face_numbereach_test_face = test_features(each_test_face_index, :)'; % 讀取測試用的臉特征% 初始化最小距離min_distance = inf;predicted_label = -1;% 對比每一張訓練臉for each_train_face_index = 1:size(train_features, 1)% 計算距離distance = norm(each_test_face - train_features(each_train_face_index, :)');% 如果當前距離小于之前的最小距離，則更新最小距離和標簽if distance < min_distancemin_distance = distance;predicted_label = floor((train_data_index(1, each_train_face_index) - 1) / j) + 1;endend% 獲取實際標簽real_label = floor((test_data_index(1, each_test_face_index) - 1) / j) + 1;% 判斷預測結果是否正確if predicted_label == real_labelcorrect_predict_number = correct_predict_number + 1;endend% 計算正確率correct_rate = correct_predict_number / test_face_number;Y(i) = correct_rate;fprintf("維度=%d，總測試樣本：%d，正確數:%d，正確率：%1f\n", d, test_face_number, correct_predict_number, correct_rate);
end
toc% 制作折線圖使正確率隨維度變化可視化
figure;
plot(10:10:40, Y, '-o');  % x軸為維度數，y軸為正確率
xlabel('維度數');
ylabel('正確率');
title('人臉識別正確率與維度變化');

B2DPCA

%1，測試與訓練圖像的錄入
tic
reshaped_faces=[];
for i=1:40for j=1:10 a=imread(strcat('C:\Users\18023\Desktop\粵海小清華\homework\計導文件\人臉數據庫\ORL56_46\orl',num2str(i),'_',num2str(j),'.bmp'));     b = reshape(a,2576,1); b=double(b);     reshaped_faces=[reshaped_faces, b];  end
end%用每個人前六張作為訓練圖，后四張作為測試圖
%定義測試集集矩陣
test_data_index = [];
%定義訓練集矩陣
train_data_index = [];  
for k=0:(size(reshaped_faces,2)/j)-1train_data_index= [train_data_index j*k+1:j*k+6];%建立一個用于查找所需訓練圖像的目錄test_data_index= [test_data_index j*k+7:j*(k+1)];%同理，是測試圖像的目錄
end
train_data = reshaped_faces(:,train_data_index);%讀取矩陣 reshaped_faces的數據添加到矩陣train_data中
test_data = reshaped_faces(:, test_data_index);%讀取矩陣 reshaped_faces的數據添加到矩陣test_data中
mean_face = mean(train_data,2);%2，協方差矩陣求取
%對于訓練矩陣重新還原成二維疊加的狀態
train_matrices = cell(size(train_data, 2), 1);
for i = 1:size(train_data, 2)train_matrices{i} = reshape(train_data(:, i), [56, 46]);
end%求所有訓練圖像矩陣的均值
mean_face_matrix = zeros(56, 46);
for i = 1:length(train_matrices)mean_face_matrix = mean_face_matrix + train_matrices{i};
end
mean_face_matrix = mean_face_matrix / length(train_matrices);% 行和列方向特征值分解
%行方向
cov_2dpca_col = zeros(46, 46);
for i = 1:length(train_matrices)centered_matrix = train_matrices{i} - mean_face_matrix;%求得每一個二維數據去中心化的結果cov_contribution = centered_matrix' * centered_matrix;%行方向協方差矩陣求解，用轉置乘以原矩陣cov_2dpca_col = cov_2dpca_col + cov_contribution;
end
cov_2dpca_col = cov_2dpca_col / length(train_matrices);
%列方向
cov_2dpca_row = zeros(56, 56);
for i = 1:length(train_matrices)centered_matrix = train_matrices{i} - mean_face_matrix;%求得每一個二維數據去中心化的結果cov_contribution = centered_matrix * centered_matrix';%列方向協方差矩陣求解，用原矩陣乘以其轉置cov_2dpca_row = cov_2dpca_row + cov_contribution;
end
cov_2dpca_row = cov_2dpca_row / length(train_matrices);% 列方向特征向量求取
[eigen_vectors_col, eigenvalues_col] = eig(cov_2dpca_col);
eigen_values_col = diag(eigenvalues_col);
[sorted_eigen_values_col, index_col] = sort(eigen_values_col, 'descend');
sorted_eigen_vectors_col = eigen_vectors_col(:, index_col);% 行方向特征向量求取
[eigen_vectors_row, eigenvalues_row] = eig(cov_2dpca_row);
eigen_values_row = diag(eigenvalues_row);
[sorted_eigen_values_row, index_row] = sort(eigen_values_row, 'descend');
sorted_eigen_vectors_row = eigen_vectors_row(:, index_row);%3，人臉重構查看效果
figure;
for idx = 1:8dim_col = idx * 5;  % 列方向降維dim_row = idx * 5;  % 行方向降維% 行和列方向特征向量選取eigen_faces_col = sorted_eigen_vectors_col(:, 1:dim_col);eigen_faces_row = sorted_eigen_vectors_row(:, 1:dim_row);% 獲取第一張訓練圖像并重構test_image = train_matrices{1};centered_image = test_image - mean_face_matrix;% B2DPCA投影和重構projection = eigen_faces_row' * centered_image * eigen_faces_col;rebuild_face = mean_face_matrix + eigen_faces_row * projection * eigen_faces_col';% 顯示重構結果subplot(2, 4, idx);imshow(mat2gray(rebuild_face));title(sprintf('維度=%d×%d', dim_row, dim_col));
end%4,人臉識別部分
index = 1;       
k = length(10:10:40);
Y = zeros(1, k);% 將測試圖像轉換為矩陣形式
test_matrices = cell(size(test_data, 2), 1);
for i = 1:size(test_data, 2)test_matrices{i} = reshape(test_data(:, i), [56, 46]);
end% 遍歷每個維度
for i = 1:kd = 10 * i;  % 每個方向的維度% 取出特征向量projection_vectors_col = sorted_eigen_vectors_col(:, 1:d);projection_vectors_row = sorted_eigen_vectors_row(:, 1:d);% 計算訓練和測試圖像的特征train_features = zeros(size(train_data, 2), d * d);test_features = zeros(size(test_data, 2), d * d);% 對每張訓練圖像進行雙向投影for j = 1:length(train_matrices)centered_matrix = train_matrices{j} - mean_face_matrix;projected = projection_vectors_row' * centered_matrix * projection_vectors_col;train_features(j, :) = projected(:)';end% 對每張測試圖像進行雙向投影for j = 1:length(test_matrices)centered_matrix = test_matrices{j} - mean_face_matrix;projected = projection_vectors_row' * centered_matrix * projection_vectors_col;test_features(j, :) = projected(:)';end% 測試前準備過程correct_predict_number = 0;test_face_number = size(test_features, 1);% 遍歷每一個待測試人臉for each_test_face_index = 1:test_face_numbereach_test_face = test_features(each_test_face_index, :)'; % 讀取測試用的臉特征% 初始化最小距離min_distance = inf;predicted_label = -1;% 對比每一張訓練臉for each_train_face_index = 1:size(train_features, 1)% 計算距離distance = norm(each_test_face - train_features(each_train_face_index, :)');% 如果當前距離小于之前的最小距離，則更新最小距離和標簽if distance < min_distancemin_distance = distance;predicted_label = floor((train_data_index(1, each_train_face_index) - 1) / j) + 1;endend% 獲取實際標簽real_label = floor((test_data_index(1, each_test_face_index) - 1) / j) + 1;% 判斷預測結果是否正確if predicted_label == real_labelcorrect_predict_number = correct_predict_number + 1;endend% 計算正確率correct_rate = correct_predict_number / test_face_number;Y(i) = correct_rate;fprintf("維度=%d，總測試樣本：%d，正確數:%d，正確率：%1f\n", d, test_face_number, correct_predict_number, correct_rate);
end
toc% 制作折線圖使正確率隨維度變化可視化
figure;
plot(10:10:40, Y, '-o');  % x軸為維度數，y軸為正確率
xlabel('維度數');
ylabel('正確率');
title('人臉識別正確率與維度變化');