matlab的Feature Transformation提供了一个有关主成分分析的介绍和例子。php
多元统计分析中广泛存在的困难中,有一个困难是多元数据的可视化。matlab的plot能够显示两个变量之间的关系,plot3和surf能够显示三维的不一样。可是当有多于3个变量时,要可视化变量之间的关系就很困难了。函数
幸运的是,在一组多变量的数据中,不少变量经常是一块儿变更的。一个缘由是不少变量是同一个驱动影响的的结果。在不少系统中,只有少数几个这样的驱动,可是多余的仪器使咱们测量了不少的系统变量。当这种状况发生的时候,你须要处理的就是冗余的信息。而你能够经过用一个简单的新变量代替这组变量来简化此问题。ui
主成分分析是一个定量的严格的能够起到简化做用的方法。它产生一组叫作主成分的新变量,每个主成分是原始变量的线性组合。全部主成分是相互正交的,从而不存在冗余的信息。全部主成分造成了原始数据空间的一组正交基。spa
可是有无数种方式来建立一组正交基,那主成分正交基有什么特殊之处呢?3d
第一主成分是数据空间的一个轴。当你把各观察值投影到这个轴时,结果会造成一个新变量,这个变量的方差是全部可选的轴中最大的。orm
第二主成分是空间的另外一个轴,它垂直于第一个轴。投影数据到这个轴上将获得另外一组变量,此变量的方差是全部可选的轴中最大的。blog
最后获得的全部主成分个数是与原始变量相同的,可是一般前几个主成分方差的和占到了原始数据总方差的80%以上。经过绘制这组新变量,研究者经常会更深刻的了解产生原始数据的驱动力。ip
在matlab中,可使用princomp函数来找到主成分,经过使用此函数,你须要用到测量的原始数据。然而,若是你缺乏这些数据,可是有这些数据的样本自相关或协方差矩阵,你能够经过使用pcacov函数来作主成分分析。ci
考虑一个例子,它使用了衡量美国329个城市生活质量的9个指标:气候、住房、健康、犯罪率、交通、教育、艺术、娱乐和经济。对于各指标,越高表示越好,如高的犯罪指标表示低的犯罪率。rem
先从cities.mat中加载数据
>> load cities
>> whos
Name Size Bytes Class Attributes
categories 9x14 252 char
names 329x43 28294 char
ratings 329x9 23688 double
cities.mat中包含了3个变量:
categories--一个包含各指标名字的字符串矩阵
names--一个包含329个城市名字的字符串矩阵
ratings--329行9列的数据矩阵
为快速对ratings数据有个直观的印象,可使用盒图绘制出数据来
boxplot(ratings,'orientation','horizontal','labels',categories)
能够注意到大致上艺术和住房的指标要比犯罪和睦候的指标有更大的变异性。
一般你会考虑绘制各对原始变量,但这将有36对变量图。而主成分分析或许会下降你须要考虑的变量个数。
有时对原始数据计算主成分是能够的,这在各变量是相同的单位时是很合适的。可是当变量单位不一样或不一样列数据的方差相差很大时,先对数据作标准化是更好的。
能够经过除以各列标准差来标准化数据,以下
stdr = std(ratings);
sr = ratings./repmat(stdr,329,1);
这样就为找主成分作好准备了。
[coefs,scores,variances,t2] = princomp(sr);
下面各部分解释了princomp函数的输出。
coefs包含了产生主成分的原始变量线性组合的系数,它常被称为loadings。它是9X9的矩阵,前三个主成分的系数向量为
>> c3 = coefs(:,1:3)
c3 =
0.2064 -0.2178 0.6900
0.3565 -0.2506 0.2082
0.4602 0.2995 0.0073
0.2813 -0.3553 -0.1851
0.3512 0.1796 -0.1464
0.2753 0.4834 -0.2297
0.4631 0.1948 0.0265
0.3279 -0.3845 0.0509
0.1354 -0.4713 -0.6073
能够看到在第一主成分中最大的系数是第三个和第七个,即健康和艺术。此主成分的全部系数都是正的。
主成分一般是单位长且正交的
I = c3'*c3
I =
1.0000 -0.0000 -0.0000
-0.0000 1.0000 -0.0000
-0.0000 -0.0000 1.0000
scores包含了原始数据在被主成分定义的新坐标系统中的坐标,它是和原始数据矩阵同大小的。
绘制scores的前两列能够展现原始数据投影到前两个主成分的新数据。princomp函数计算的scores具备零均值。
plot(scores(:,1),scores(:,2),'+')
xlabel('1st Principal Component')
ylabel('2nd Principal Component')
注意到上图中右边的离群点。函数gname能够用于标出这些离群点,调用gname时传入一个string矩阵,它包含了各点的标签,以下
gname(names)
移动鼠标在右半部分的点点击,当你点击各点的时候,在点上会标记names中对应的字符串,以下
从上面的标记,发现这些离散点是美国一些人口比较多的城市,它们明显与其余数据不一样,因此它们可能须要被分离开。为移除这些数据,首先识别这些数据的行,方法以下
1.关闭上面的figure
2.重绘plot
plot(scores(:,1),scores(:,2),'+')
xlabel('1st Principal Component');
ylabel('2nd Principal Component');
3.运行gname函数,如输入参数
4.标记离散点,标记自动为这些数据的行数
而后建立一个包含这些点的变量
metro = [43 65 179 213 234 270 314];
names(metro,:)
ans =
Boston, MA
Chicago, IL
Los Angeles, Long Beach, CA
New York, NY
Philadelphia, PA-NJ
San Francisco, CA
Washington, DC-MD-VA
而后移除这些行
rsubset = ratings;
nsubset = names;
nsubset(metro,:) = [];
rsubset(metro,:) = [];
size(rsubset)
ans =
322 9
variances是一个包含了主成分方差的向量,scores的每一列的方差即为variances对应的项
variances
variances =
3.4083
1.2140
1.1415
0.9209
0.7533
0.6306
0.4930
0.3180
0.1204
你能够很容易的计算全部差别性占到百分比
percent_explained = 100*variances/sum(variances)
percent_explained =
37.8699
13.4886
12.6831
10.2324
8.3698
7.0062
5.4783
3.5338
1.3378
使用pareto函数能够以图的方式显示此百分比
pareto(percent_explained)
xlabel('Principal Component')
ylabel('Variance Explained (%)')
显示的figure中能够看到第一主成分明显比第二主成分要高不少,但它只解释了不到40%的方差,因此须要使用更多的成分。能够看到前三个主成分解释了三分之二的方差,因此这三个成分能够做为一个理想的方式来下降维数。
t2是一个衡量各观测值离数据中心距离的统计量,这是一个找出数据中极值点的有效分析方法。
[st2, index] = sort(t2,'descend'); % Sort in descending order.
extreme = index(1)
extreme =
213
names(extreme,:)
ans =
New York, NY
这个极值点是不奇怪的,由于Vew York的ratings是离美国城镇平均值是最远的。
使用biplot函数能够帮助咱们可视化各变量的主成分的系数和各观察值的主成分scores。以下,
biplot(coefs(:,1:2), 'scores',scores(:,1:2),...
'varlabels',categories);
axis([-.26 1 -.51 .51]);
上图中,9个变量以向量的形式呈现出来,向量的方向和长度表示了各变量对两个主成分的贡献。例如,对于第一主成分,对全部9个变量的系数都是正值。对于第二主成分,变量教育、健康、艺术和交通是正的贡献,其余5个是负的贡献。这代表此成分在不一样城市间是有区别的:大的值的变量、小的值的变量和负值的变量。
上图中变量标签有时候会很拥挤,能够在绘图时在varlabels参数中忽略一些,或在figure的编辑模式下选中并移动它们。可使用Tools下的Data Cursors功能查看图中的信息
biplot函数也能够用于绘制三维的信息
biplot(coefs(:,1:3), 'scores',scores(:,1:3),...
'obslabels',names);
axis([-.26 1 -.51 .51 -.61 .81]);
view([30 40]);