前往小程序,Get更优阅读体验!
立即前往
首页
学习
活动
专区
工具
TVP
发布
社区首页 >专栏 >R语言 主成分分析PCA(绘图+原理)

R语言 主成分分析PCA(绘图+原理)

作者头像
拴小林
发布2021-05-31 16:21:58
13.8K0
发布2021-05-31 16:21:58
举报
文章被收录于专栏:数据驱动实践
代码语言:javascript
复制
# install.packages("devtools")
# devtools::install_github("vqv/ggbiplot")
library(ggbiplot)
data(wine)
head(wine)
# > head(wine)
# Alcohol MalicAcid  Ash AlcAsh  Mg Phenols Flav NonFlavPhenols Proa Color  Hue   OD Proline
# 1   14.23      1.71 2.43   15.6 127    2.80 3.06           0.28 2.29  5.64 1.04 3.92    1065
# 2   13.20      1.78 2.14   11.2 100    2.65 2.76           0.26 1.28  4.38 1.05 3.40    1050
# 3   13.16      2.36 2.67   18.6 101    2.80 3.24           0.30 2.81  5.68 1.03 3.17    1185
# 4   14.37      1.95 2.50   16.8 113    3.85 3.49           0.24 2.18  7.80 0.86 3.45    1480
# 5   13.24      2.59 2.87   21.0 118    2.80 2.69           0.39 1.82  4.32 1.04 2.93     735
# 6   14.20      1.76 2.45   15.2 112    3.27 3.39           0.34 1.97  6.75 1.05 2.85    1450

wine.pca <- prcomp(wine, scale. = TRUE)
## bioplot
ggbiplot(wine.pca, obs.scale = 1, var.scale = 1,
         groups = wine.class, ellipse = TRUE, circle = TRUE) +
  scale_color_discrete(name = '') +
  theme_light()

ggbiplot Arguments

代码语言:javascript
复制
ggbiplot(pcobj, choices = 1:2, scale = 1, pc.biplot =
TRUE, obs.scale = 1 - scale, var.scale = scale, groups =
NULL, ellipse = FALSE, ellipse.prob = 0.68, labels =
NULL, labels.size = 3, alpha = 1, var.axes = TRUE, circle
= FALSE, circle.prob = 0.69, varname.size = 3,
varname.adjust = 1.5, varname.abbrev = FALSE, ...)
pcobj # prcomp()或princomp()返回结果
choices # 选择轴,默认1:2
scale # covariance biplot (scale = 1), form biplot (scale = 0). When scale = 1, the inner product between the variables approximates the covariance and the distance between the points approximates the Mahalanobis distance.
obs.scale # 标准化观测值
var.scale # 标准化变异
pc.biplot # 兼容 biplot.princomp()
groups # 组信息,并按组上色
ellipse # 添加组椭圆
ellipse.prob # 置信区间
labels #向量名称
labels.size #名称大小
alpha  #点透明度 (0 = TRUEransparent, 1 = opaque)
circle #绘制相关环(only applies when prcomp was called with scale = TRUE and when var.scale = 1)
var.axes  #绘制变量线-菌相关
varname.size  #变量名大小
varname.adjust #标签与箭头距离 >= 1 means farther from the arrow
varname.abbrev # 标签是否缩写
代码语言:javascript
复制
library(psych)
library(ggplot2)
PC <- principal(wine, nfactors=2, rotate ="none")
pc <- data.frame(PC$scores)
p  <- ggplot(pc, aes(x=PC1, y=PC2,color=wine.class )) +
  geom_point(size=4,alpha=0.5)+
  theme(axis.text= element_text(size=20))+
  theme(panel.grid.major =element_blank(),
        panel.grid.minor = element_line(color="steelblue"),
        panel.background = element_blank(),
        axis.line = element_line(colour = "black"))+
  stat_ellipse(lwd=1,level = 0.8)
代码语言:javascript
复制
data(wine)
wine.pca <- prcomp(wine, scale. = TRUE)
## bioplot
TB= data.frame(wine.pca$x)
TB$class = wine.class

p0 <- ggbiplot(wine.pca, obs.scale = 1, var.scale = 1,
               groups = wine.class, ellipse = TRUE, circle = TRUE) +
  scale_color_discrete(name = '') +
  theme_light()+theme(legend.position = 'none',
                      axis.title=element_blank())

p1 <- ggplot() +
  geom_density(data=TB, aes(x=PC1,fill=class),alpha=0.5) +
  theme_map() +
  theme(legend.position = 'none')

p2 <- ggplot() +
  geom_density(data=TB, aes(x=PC2,fill=class),alpha=0.5)  +
  theme_map() +theme(legend.position = 'none') +
  coord_flip()  +  scale_y_reverse()

p3 <- ggplot(TB, aes(x=class,y=PC1)) +
  geom_tufteboxplot(aes(color=class),median.type = "line",
                    hoffset = 0, width = 3) +
  coord_flip()+ theme_map()+
  theme(legend.position = 'none')

p4 <- ggplot(TB, aes(x=class,y=PC2)) +
  geom_tufteboxplot(aes(color=class),median.type = "line",
                    hoffset = 0, width = 3) +
  theme_map()+
  theme(legend.position = 'none')


TB= data.frame(wine.pca$x)
TB$class = wine.class

GGlay <-
  "#BBBBBB#
CAAAAAAE
CAAAAAAE
CAAAAAAE
CAAAAAAE
CAAAAAAE
CAAAAAAE
#DDDDDD#"
p0+p1+p2+p3+p4+plot_layout(design = GGlay)

PCA原理

PCA 是一种较为常用的降维技术,PCA 的思想是将n维特征映射到k维上,这k维是全新的正交特征。这k维特征称为主元,是重新构造出来的k维特征。在 PCA 中,数据从原来的坐标系转换到新的坐标系下,新的坐标系的选择与数据本身是密切相关的。其中,第一个新坐标轴选择的是原始数据中方差最大的方向,第二个新坐标轴选取的是与第一个坐标轴正交且具有最大方差的方向,依次类推,我们可以取到这样的k个坐标轴。

PCA 操作流程

  1. 去均值,即每一位特征减去各自的平均值(当然,为避免量纲以及数据数量级差异带来的影响,先标准化是必要的)
  2. 计算协方差矩阵
  3. 计算协方差矩阵的特征值与特征向量
  4. 对特征值从大到小排序
  5. 保留最大的k个特征向量
  6. 将数据转换到k个特征向量构建的新空间中

1. 常用术语

(1)标准化(Scale)

如果不对数据进行scale处理,本身数值大的基因对主成分的贡献会大。如果关注的是变量的相对大小对样品分类的贡献,则应SCALE,以防数值高的变量导入的大方差引入的偏见。但是定标(scale)可能会有一些负面效果,因为定标后变量之间的权重就是变得相同。如果我们的变量中有噪音的话,我们就在无形中把噪音和信息的权重变得相同,但PCA本身无法区分信号和噪音。在这样的情形下,我们就不必做定标。

(2)特征值 (eigen value)

特征值与特征向量均为矩阵分解的结果。特征值表示标量部分,一般为某个主成分的方差,其相对比例可理解为方差解释度或贡献度 ;特征值从第一主成分会逐渐减小。

(3)特征向量(eigen vector)

特征向量为对应主成分的线性转换向量(线性回归系数),特征向量与原始矩阵的矩阵积为主成分得分。特征向量是单位向量,其平方和为1。特征向量主要起转换作用,其数值不能说明什么问题,解释力更强的是载荷loadings,但很多R输出中经常混用,engen vector与loadings。

(4)载荷(loading)

因子载荷矩阵并不是主成分的特征向量,即不是主成分的系数。主成分系数的求法:各自因子载荷向量除以各自因子特征值的算数平方根。特征向量是单位向量,特征向量乘以特征值的平方根构造了载荷loading。列上看,不同变量对某一PC的loadings的平方和等于其征值,因此每个变量的loadings值可表征其对PC的贡献。行上看,同一变量对不同PCs的loadings行平方和为1,表征不同PCs对某一变量方差的解释度。

(5)得分(score)

指主成分得分,矩阵与特征向量的积。·

2. PCA分析过程

手动计算

代码语言:javascript
复制
library(dplyr)
#特征分解
dat_eigen<-scale(iris[,-5],scale=T)%>%cor()%>%eigen()
#特征值提取
dat_eigen$values 
#特征向量提取
dat_eigen$vectors
#求loadings=eigen vector*sqrt(eigen value),与princomp不同
#主成分载荷表示各个主成分与原始变量的相关系数。
sweep(dat_eigen$vectors,2,sqrt(dat_eigen$values),"*")
#将中心化的变量矩阵得到每个观测值的得分
scale(iris[,-5],scale=T)%*%dat_eigen$vectors%>%head()

2.1 prcomp函数

prcomp函数使用较为简单,但是不同于常规的求取特征值和特征向量的方法,prcomp函数是对变量矩阵(相关矩阵)采用SVD方法计算其奇异值(原理上是特征值的平方根),函数帮助中描述为函数结果中的sdev。prcomp函数输入参数为变量矩阵(x),中心化(center,默认为true),标准化(scale,默认为false,建议改为true),主成份个数(rank)。prcomp函数输出有sdev(各主成份的奇异值),rotation(特征向量,回归系数),x(score得分矩阵)。

代码语言:javascript
复制
iris.pca<-prcomp(iris[,-5],scale=T,rank=4,retx=T) #相关矩阵分解
#retx表四返回score,scale表示要标准化
summary(iris.pca) #方差解释度
iris.pca$sdev #特征值的开方
iris.pca$rotation #特征向量,回归系数
iris.pca$x #样本得分score

2.2 princomp函数

princomp以计算相关矩阵或者协方差矩阵的特征值为主要手段。princomp函数输出有主成份的sd,loading,score,center,scale.prcomp函数使用较为简单,但是不同于常规的求取特征值和特征向量的方法,prcomp函数是对变量矩阵(相关矩阵)采用SVD方法计算其奇异值(原理上是特征值的平方根),函数帮助中描述为函数结果中的sdev。prcomp函数输入参数为变量矩阵(x),中心化(center,默认为true),标准化(scale,默认为false,建议改为true),主成份个数(rank)。prcomp函数输出有sdev(各主成份的奇异值及其方差累积),rotation(载荷矩阵),x(得分矩阵),center(变量的均值),scale(变量的标准偏差)

代码语言:javascript
复制
data(wine) #三种葡萄酿造的红酒品质分析数据集
wine.pca<-princomp(wine,cor=T,scores=T) 
#默认方差矩阵(cor=F),改为cor=T则结果与prcomp相同
summary(wine.pca) #各主成份的SVD值以及相对方差
wine.pca$loading #特征向量,回归系数
wine.pca$score
screenplot(wine.pca) #方差分布图
biplot(wine.pca,scale=F) #碎石图,直接把x与rotation绘图,而不标准化

2.3 psych::principal

实际上该principal主要用于因子分析。

代码语言:javascript
复制
model_pca<-psych::principal(iris[,-5],nfactors=4,rotate="none")
model_pca$values # 特征值=sdev^2
# 此处loadings与手动计算相同=prcomp的rotation*sdev 
model_pca%>%.$loadings #载荷,不是特征向量
#此处score=prcomp的score/sdev
model_pca$scores[1:5,] #此处为因子得分,不是主成分得分
model_pca$weights #此处为上面loadings/特征值,也称成份得分系数或者因子系数

3. PCA结果解释

下文引用chentong的内容

代码语言:javascript
复制
prcomp函数会返回主成分的标准差、特征向量和主成分构成的新矩阵。
不同主成分对数据差异的贡献和主成分与原始变量的关系。
1. 主成分的平方为为特征值,其含义为每个主成分可以解释的数据差异,计算方式为 eigenvalues = (pca$sdev)^2
2. 每个主成分可以解释的数据差异的比例为 percent_var = eigenvalues*100/sum(eigenvalues)
3. 可以使用summary(pca)获取以上两条信息。

这两个信息可以判断主成分分析的质量:
    成功的降维需要保证在前几个为数不多的主成分对数据差异的解释可以达到80-90%。

指导选择主成分的数目:
  1.  选择的主成分足以解释的总方差大于80% (方差比例碎石图)
  2. 从前面的协方差矩阵可以看到,自动定标(scale)的变量的方差为1 (协方差矩阵对角线的值)。待选择的主成分应该是那些方差大于1的主成分,即其解释的方差大于原始变量(特征值碎石图,方差大于1,特征值也会大于1,反之亦然)。

鉴定核心变量和变量间的隐性关系:
    原始变量与主成分的相关性 Variable correlation with PCs (var.cor) = loadings * sdev
    原始数据对主成分的贡献度 var.cor^2 / (total var.cor^2)

原文链接:https://blog.csdn.net/nikang3148/article/details/85246555

本文参与 腾讯云自媒体同步曝光计划,分享自微信公众号。
原始发表:2021-05-08,如有侵权请联系 cloudcommunity@tencent.com 删除

本文分享自 数据驱动实践 微信公众号,前往查看

如有侵权,请联系 cloudcommunity@tencent.com 删除。

本文参与 腾讯云自媒体同步曝光计划  ,欢迎热爱写作的你一起参与!

评论
登录后参与评论
0 条评论
热度
最新
推荐阅读
目录
  • 3. PCA结果解释
领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档