数学建模笔记-第十四讲-主成分分析

主成分分析

经典的降维算法PCA

数据降维

主成分分析思想

PCA计算过程

• 标准化：去除量纲
• 归一化：仅是数据映射到[0,1]上，在评价问题中方便解释

• 按列进行标准化

• 对标准化后的矩阵计算协方差矩阵

  • 第五章学过，标准化之后求协方差，消除了量纲的影响，就相当于原来样本的相关系数

• 当然也可以两步合成一步，直接计算x矩阵的样本相关系数矩阵（在第五章皮尔逊相关系数那讲提过）

• R是半正定矩阵有严格的数学推导
• 特征值求和等于矩阵的迹，这里就是p

花了10分钟，终于弄懂了特征值和特征向量到底有什么意义

深入理解矩阵特征值与特征向量的物理意义

学完线性代数一年多了都快忘记啥是特征值和特征向量了

截取：

求特征向量，就是把矩阵A所代表的空间进行正交分解，使得A的向量集合可以表示为每个向量a在各个特征向量上的投影长度。我们通常求特征值和特征向量即为求出这个矩阵能使哪些向量只发生拉伸，而方向不发生变化，观察其发生拉伸的程度。这样做的意义在于，看清一个矩阵在哪些方面能产生最大的分散度（scatter），减少重叠，意味着更多的信息被保留下来

每一个贡献率就是特征值/特征值的和

算累计贡献率

确定主成分：比如前两个主成分累积贡献率超过80%，那就取前两个主成分

主成分的系数是其特征向量

• 用于聚类
  • 本来有十个指标，无法展示
  • 主成分分析成两个指标的话，就可以
• 回归分析
  • 降维，来缓解多重共线性问题

主成分分析的应用

例1

• 这里数据是样本相关系数矩阵
• 那么相当于直接就到了第二步

• matlab计算的特征向量默认已经归一化，模长是1。即平方和开根号为1
• 特征值和为指标个数：6
• 贡献率=特征值/指标个数
• 根据累积贡献率，考虑前面三个成分

• X时已经标准化的指标
• 载荷就可以理解为这里指标前面的系数
• 根据分析，可以判断三个主成分的特点。

主成分的说明

• 困难就是要对主成分给出较好的解释
• 要取舍（如例1舍弃第三个主成分），还有能否解释运气

例2

这里其实指标也不算多，就8个，还是可以写全的。如果说太多指标的话，可以像上面一样不写全

• 解释F1和F2 （前两个主成分）

• F1:田径上有强有弱，反应了强弱成分
• F2:反应速度和耐力成绩的对比

MATLAB

主成分分析是对指标进行线性组合，适合对相关性强的指标进行降维

• 得到样本相关系数矩阵

clear;clc
 load data1.mat   % 主成分聚类
%  load data2.mat   % 主成分回归

% 注意，这里可以对数据先进行描述性统计
% 描述性统计的内容见第5讲.相关系数
[n,p] = size(x);  % n是样本个数，p是指标个数

%% 第一步：对数据x标准化为X
X=zscore(x);   % matlab内置的标准化函数（x-mean(x)）/std(x)

%% 第二步：计算样本协方差矩阵
R = cov(X);

%% 注意：以上两步可合并为下面一步：直接计算样本相关系数矩阵
R = corrcoef(x);
disp('样本相关系数矩阵为：')
disp(R)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18

• 计算R的特征值和特征向量

%% 第三步：计算R的特征值和特征向量
% 注意：R是半正定矩阵，所以其特征值不为负数
% R同时是对称矩阵，Matlab计算对称矩阵时，会将特征值按照从小到大排列哦
% eig函数的详解见第一讲层次分析法的视频
[V,D] = eig(R);  % V 特征向量矩阵  D 特征值构成的对角矩阵

%  rot90函数可以使一个矩阵逆时针旋转90度，然后再转置，就可以实现将矩阵的列颠倒的效果
V=rot90(V)';
1
2
3
4
5
6
7
8

• diag()函数

  • 返回主对角线上的元素

• 计算主成分贡献率和累积贡献率

  • 由于我们要把特征值从大到小排序，使用lambda = lambda(end:-1:1);

  • 使特征向量矩阵也要跟特征值对应（最后一列变为第一列，倒数第二列变为第二列等），进行颠倒，比较巧妙：先逆时针旋转90度，再转置

  • %  rot90函数可以使一个矩阵逆时针旋转90度，然后再转置，就可以实现将矩阵的列颠倒的效果
    V=rot90(V)';
    1
    2

  • cumsum()是求累加值的函数

%% 第四步：计算主成分贡献率和累计贡献率
lambda = diag(D);  % diag函数用于得到一个矩阵的主对角线元素值(返回的是列向量)
lambda = lambda(end:-1:1);  % 因为lambda向量是从小大到排序的，我们将其调个头
contribution_rate = lambda / sum(lambda);  % 计算贡献率
cum_contribution_rate = cumsum(lambda)/ sum(lambda);   % 计算累计贡献率  cumsum是求累加值的函数
disp('特征值为：')
disp(lambda')  % 转置为行向量，方便展示
disp('贡献率为：')
disp(contribution_rate')
disp('累计贡献率为：')
disp(cum_contribution_rate')
disp('与特征值对应的特征向量矩阵为：')
% 注意：这里的特征向量要和特征值一一对应，之前特征值相当于颠倒过来了，因此特征向量的各列需要颠倒过来
%  rot90函数可以使一个矩阵逆时针旋转90度，然后再转置，就可以实现将矩阵的列颠倒的效果
V=rot90(V)';
disp(V)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17

• 计算我们需要的主成分的值
  • 也就是每个样本每个指标的值**(x)，乘上该主成分特征向量的对应系数(a)**

%% 计算我们所需要的主成分的值
m =input('请输入需要保存的主成分的个数:  ');
F = zeros(n,m);  %初始化保存主成分的矩阵（每一列是一个主成分）
for i = 1:m
    ai = V(:,i)';   % 将第i个特征向量取出，并转置为行向量
    Ai = repmat(ai,n,1);   % 将这个行向量重复n次，构成一个n*p的矩阵
    F(:, i) = sum(Ai .* X, 2);  % 注意，对标准化的数据求了权重后要计算每一行的和
end
1
2
3
4
5
6
7
8

对结果的解释

• 综合消费支出
• 观察到衣着和医疗的正载荷较高，消费倾向成分（有一些经济学知识）

这里如果给不出合理的解释，对后续的分析较难进行，所以说主成分分析最难的步骤其实就在于此

主成分分析的滥用：主成分得分

• 像topsis那样不会损失信息，已经有数据了没必要用主成分分析

• 把指标的属性改了

主成分分析用于聚类

计算出第一主成分和第二主成分的值，将其视为两个新的指标
（可以在图上直观的展示各样本的分布情况）

将得到的F1，F2导入spss聚类分析

%%主成分聚类 ： 将主成分指标所在的F矩阵复制到Excel表格，然后再用Spss聚类
% 在Excel第一行输入指标名称（F1,F2, ..., Fm）
% 双击Matlab工作区的F,进入变量编辑中，然后复制里面的数据到Excel表格
% 导出数据之后，我们后续的分析就可以在Spss中进行。
1
2
3
4
5

系统聚类，勾选上谱系图

分为三类：广东、（上海浙江北京）、其他

再系统聚类设定为三类，得到分析结果，更具结果画图（具体操作看第十讲）

• 跟第十讲的聚类结果是有一些差别的，因为主成分分析会损失信息，但最大的意义是可视化
• 根据图可以大概看出来聚类的效果
• 只有指标个数特别多，指标之间相关性很强（第一第二主成分的累积贡献率就会比较大，损失的就会少）

主成分回归

• 例子

• 当时逐步回归，可能只用了三个指标，其他的没用上。这里用主成分分析对指标进行线性组合，虽然也会有损失，但比没用损失少。
• 更换数据运行matlab代码

• 论文中解释一下前两个主成分，要好好琢磨
• 进行异方差检验
• 用F1、F2进行回归，解释回归系数。
  • 因变量y也要标准化
  • 综合性投入每增加1个单位，对产量会增加。。

说明

• 如果能解释清楚主成分，那么可以既用主成分分析，也用逐步回归
• 如果解释不清，还是用逐步回归

作业

本笔记来自清风老师的数学建模，强烈推荐该课程！