试使用matlab的PCA函数对Yale人脸数据进行降维，并观察前20个特征向量对应的图像

学习笔记：
PCA(Principal Component Analysis)，即主成分分析方法，是一种使用最广泛的数据降维算法。PCA的主要思想是将n维特征映射到k维上，
这k维是全新的正交特征也被称为主成分，是在原有n维特征的基础上重新构造出来的k维特征。PCA的工作就是从原始的空间中顺序地找一组相互正交的坐标轴，新的坐标轴的选择与数据本身是密切相关的。其中：
第一个新坐标轴选择是原始数据中方差最大的方向，
第二个新坐标轴选取是与第一个坐标轴正交的平面中使得方差最大的，
第三个轴是与第1,2个轴正交的平面中方差最大的。
依次类推，可以得到n个这样的坐标轴。通过这种方式获得的新的坐标轴，我们发现，大部分方差都包含在前面k个坐标轴中，后面的坐标轴所含的方差几乎为0。于是，我们可以忽略余下的坐标轴，只保留前面k个含有绝大部分方差的坐标轴。事实上，这相当于只保留包含绝大部分方差的维度特征，而忽略包含方差几乎为0的特征维度，实现对数据特征的降维处理。

matlab直接调用pca即可，关键是pca函数的输入输出参数的理解。

第一种理解方式：
[coeff,score,latent,tsquare]=pca(X)

首先，输入X是一个 n×p 的矩阵，在这里，n是样本个数，p是原始的维数。针对这里的图像数据集，我们可以认为每一幅图像对应一个维度，
而将二维图像展开成一维向量后，每个元素位置对应一个样本。比如图像大小为 243*320，一共有165张。则样本数就是77760，维数就是165。
所以输入的X的大小就是 77760×165 的形状。

至于这里为什么将每一幅图像认为是一个维度，个人理解是，我们的目标是用较少的脸可以尽可能准确地组合重建出原始的165个脸，也就是说要降低脸的个数，
也就是要降维的意思。所以这里认为每一个脸为一个维度，然后进行降维。而将图像转成向量主要是因为输入限制。

coeff,解释为主成分系数，其实就是输入的 X 的协方差矩阵 对应的所有特征向量组成的矩阵，并且根据每个特征向量对应的特征值，按列递减的顺序排列。

score是对主成分的打分，也就是原X矩阵在主成分空间的表示，每行对应一个样本观测值，每一列对应一个主成分，行列数目和输入的X的一致。

latent是一个向量，就是 X 的协方差矩阵对应的特征值组成的向量，按递减顺序排列。

tsquare是表示对每个样本点的Hotelling的T方统计量。

结果分析：
database为原始灰度图，result为PCA降维（k=20）后的灰度图。
可以看出特征少时有很大的重影，但由于是最显著的几种特征，依然能看出与原图几分相似。
降维前图像：
降维后图像（K=20）：
降维后图像（K=100）：

降维后图像（k=150）：

此外，可以发现当k从20增加到150时，选择的k个特征的累计方差占比已经接近于1了，
而相应地，图像的特征已经与原始图像非常接近了，这意味着我们可以用大约150维的向量来描述一张原本纬数达70000多维的图像，
可见PCA在这样的灰度人脸图像下的降维是非常有效的。

代码清单：

name_prefix1 = 'database//s%d.bmp';  %加载未处理图像
name_prefix2 = 'result//s%d.bmp';    %保存PCA降维后图像
[k, n] = deal(20, 165);              %取前20个特征向量，共165张图片
M = read_image(name_prefix1, n);
M = dimension_reduction(M, k, n);
write_image(name_prefix2, M, n)
 
function [M] = read_image(name_prefix, n)
    %读取
    M = zeros(100*100, n);
    for i=1:n 
        file_name = sprintf(name_prefix, i);
        image = imread(file_name, 'bmp');
        imshow(image);      % 显示降维前的灰度图
        M(:,i) = image(:);  % M矩阵的每一行为一张图像矩阵
    end
end
function [M] = dimension_reduction(M, k, n)
    %PCA分析
    F = pca(M);        % F为协方差矩阵的所有特征向量降序排列组成的矩阵
    F(:, k+1:n) = 0;   % 只保留前k个特征
    %F = F(:,1:k);
    M = M * (F * F');  % F*F'是散度矩阵
end
function [] = write_image(name_prefix, M, n)
    %写入
    for i=1:n
        file_name = sprintf(name_prefix, i);
        img = M(:,i);
        maxx = max(img);
        minn = min(img);
        img = 255*(img-minn)/(maxx-minn);     % 将值转到0~255
        img = reshape(M(:,i), 100, 100);
        img = mat2gray(img);      % 从数据矩阵创建灰度图
        imshow(img);              % 显示降维后的灰度图
        imwrite(img, file_name, 'bmp');
    end
end