数据降维
在进行数据挖掘或者机器学习时,我们面临的数据往往是高维数据。
相较于低维数据,高维数据为我们提供了更多的信息和细节,也更好的描述了样本;但同时,很多高效且准确的分析方法也将无法使用。
处理高维数据和高维数据可视化是数据科学家们必不可少的技能。解决这个问题的方法便是降低数据的维度。在数据降维时,要使用尽量少的维度来表达较多原数据的特性和结构。
降维是通过减少数据中的指标(或变量)以化简数据的过程。这里的减少指标,并不是随意加减,而是用复杂的数理知识,得到几个“综合指标”来代表整个数据。
PCA原理
主成分分析(Principal Component Analysis, PCA)是一种线性降维算法,也是一种常用的数据预处理(Pre-Processing)方法。
它的目标是是用方差(Variance)来衡量数据的差异性,并将差异性较大的高维数据投影到低维空间中进行表示。
绝大多数情况下,我们希望获得两个主成分因子:分别是从数据差异性最大和次大的方向提取出来的,称为PC1(Principal Component 1) 和 PC2(Principal Component 2)。
PCA原理简述
PCA具体是如何实现的?PC1和PC2又是如何计算的呢?
举个简单的例子,Scores.xlsx 包含了约70名学生的全科考试成绩。其中每名学生是一个独立的样本,每门学科的成绩都是一个数据维度(共有13门成绩)。我们的目的是通过分析学生的考试成绩来判断学生的类别(理科、文科生,和体育、艺术特长生)。
下表为随机选取的六名学生的语文和数学考试成绩:
Step 1 散点图
制作为散点图:
上图中每个点代表了一个学生,X轴代表语文成绩,Y轴代表数学成绩。
Step 2 计算均值
然后分别取所有样本的X平均值和Y平均值,并将这两个值变为X、Y坐标,在图中画出这个点(用五角星表示):
Step 3 中心化
按照图中箭头所示方向,将整个坐标系平移,使原点与五角星重叠。这样就获得了一个新的平面直角坐标系:
尽管此时坐标系和每个点的值都发生了变化,点与点之间的相对位置仍保持一致。
这个过程称为“中心化”。“中心化”处理的原因是,这些数字后续会参与统计运算,比如求样本方差。
Step 4 找到最优拟合线
找到这些点的最优拟合线(Line of Best Fit),也就找到了PC1,再通过原点做PC1的垂线,就找到了PC2:
处理三维数组时便会产生第三个因子(PC3),以此类推,数据的维度越大,因子的数量也就越多。当维度大于等于4的时候,我们是无法想象出图像的,但PC4确实存在;假设有x个维度,便可以做x-1条垂线,就能得到PCx。接下来要做的便是选取最能代表数据差异性的两个因子,作为PC1和PC2。
Step 5 投影(Projection)
按照下图所示,将点A投影到PC1上(六角星的位置),并计算其与原点之间的距离称为d1:
其余的五个点也做同样操作,得出d2至d5,再求这六个距离的平方和,称为PC1的特征值(Eigenvalue)。
Step 6 计算差异值
然后将PC1的特征值除以总样本数量减一(n-1),就计算出了PC1的差异值(Variation):
Step 7 碎石图(scree plot)
以此类推,并选择差异值最大的两个因子作为PC1 和 PC2。假设在某个三维数组中,获得了PC1、PC2和PC3的差异值分别为18,7,5。通过计算(18+7)/ (18+7+5) ≈ 83.3% 得到结论:PC1 和 PC2 代表了这个三维数组83.3%的差异性。在本次分析的13个因子中,PC1和PC2描述了整组数据约81%的差异性:
Step 8 PCA图(PCA plot)
最后,再通过选中的PC1和PC2将样本映射回本身所在的坐标,就可以得到降维后的图像(PCA Plot)。
上图是在经过PCA降维处理后,得到的浮点图;图中每个点都代表了一名学生,其中文科生为绿色,理科生为红色,体育生和艺术生分别为黄色和蓝色。坐标系的X轴为PC1,Y轴为PC2。从图中我们可以明显地看出同类别学生的聚类(Clustering)趋势,也证明了PCA在降维的同时,尽可能地保留了原数据的特性。
PCA的理论推导
上述原理只是PCA原理的简单化说明,实际上,主成分不是原来的指标中的任何一个,而是由所有原有指标数据线性组合而来。
PCA是将一组变量通过正交变换转变成另一组变量的分析方法,以此来实现数据降维的目的。转换后得到的这一组变量,即是我们所说的主成分。也就是把多指标的数据用少数几个综合指标(主成分)替代,还原数据最本质特征。然后根据这些主成分对样本进行聚类,代表样本的点在坐标轴上距离越远,说明样本差异越大。同时,还可以看两个或多个分组之间能否分得开,能否找到差别,也是一种质控的手段。
这些主成分会依据方差的大小进行排序,称作主成分(PC)1、主成分2、……主成分r。而每个主成分的方差在这一组变量中的总方差中所占的比例,即是该主成分的可解释差异贡献度。一般我们要求前K个主成分的累计可解释差异贡献度大于80%,k小于4,因为最多作出三维图像,即仅考察贡献度前2或者前3的主成分,经过可视化后,即得到了二维或三维PCA散点图。在二维图像中,把横纵坐标百分比加起来,这就是它们可解释差异贡献度,这个数字越大模型的效果越好,三维图像是X+Y+Z,这个数字要求大于80%。
以最简单的模型为例:如何把二维数据压缩至一维。
数据压缩的目标是:让新数据的方差尽可能地大。这个标准能使得新数据尽可能地不丢失原有数据的信息,因为方差越大,数据间的差异越大。
如下图所示:有六个点,每个点有两个特征,分别对应x轴和y轴。需要把他们压缩成一维的数据,即每个点只有一个特征。因此要寻一条直线,让所有点投影到该直线上。该直线上的刻度即为新数据的值。
首先进行中心化处理。中心化的好处在于,我们寻求的直线必定经过原点。如下图所示,只需从所有经过原点的直线中,找一条直线,使得各个数据的方差最大。
如下图所示,设P是样本的一个数据。注意到这样一个性质:由于OP的距离是固定的,因此过原点作任何一条直线,记Q是P在该直线上的投影点,都有 OQ2 + PQ2 = 定值。
显然知道 OQ (带正负号的长度)是数据 P(xp, yp) 压缩后的值。
那么如下图所示:对于 A1 到 A6 六个点,分别作直线的垂线,交直线于 B1 到 B6 六个点。直线可以看做数轴,坐标系的原点就是数轴的原点。记 B1到 B6 六个点的刻度分别为 b1,…,b6。那么 b1,…,b6 就是降维后的数据。
也就说明 b 的均值为0。那么 b 的方差为:
![](https://img2023.cnblogs.com/blog/2975222/202304/2975222-20230401210028021-1867905178.png)
我们的目标是让压缩后的数据方差最大,即让
![](https://img2023.cnblogs.com/blog/2975222/202304/2975222-20230401205737446-109115701.png)
记点
![](https://img2023.cnblogs.com/blog/2975222/202304/2975222-20230401205843510-992469378.png)
![](https://img2023.cnblogs.com/blog/2975222/202304/2975222-20230401205856740-852493095.png)
![](https://img2023.cnblogs.com/blog/2975222/202304/2975222-20230401210003178-831402670.png)
于是
![](https://img2023.cnblogs.com/blog/2975222/202304/2975222-20230401210052584-2052328857.png)
使
![](https://img2023.cnblogs.com/blog/2975222/202304/2975222-20230401210105324-930433366.png)
![](https://img2023.cnblogs.com/blog/2975222/202304/2975222-20230401210123096-1968927820.png)
这个结果很美妙,在这里和一元线性回归作比较。如下图所示,一共有九个点。我们要进行一元线性回归,所采取的策略是让平方损失最小,通俗地来说,就是让九条绿线的平方和最小。
如果是采用PCA,如下图所示,我们采取的策略是让九条红线的平方和最小。
事实上,一元线性回归也能采取“使点到直线距离的平方和最小”的策略。两种策略各有各的优势与缺陷。只不过“平方误差最小原则”已经被广泛地接受。
类似地,对于n维数据,投影到更小的k维空间也能同样地进行类比。
为什么PCA可以代替所有数据?
主成分为什么拽到可以代替所有数据?换句话说,为什么较多的变量经过数据变换之后变成了较少的几个变量呢?这不是会有大量的信息丢失掉吗?
这是因为部分变量其实是相互关联的,这就会造成数据冗余。而降维就可以帮助我们去除这些指标中重叠、多余的信息,把数据最本质和关键的信息提取出来。在我们最开始得到的一组变量中,变量之间并不是完全相互独立的。例如一个位点发生了变异,那么与之连锁的几个位点也大概率会发生变异;或者一个基因的表达量发生了变化,同一条通路中的其他基因的表达量也大概率会发生变化,即变量之间是存在相关性的。极端一点,假设两个位点完全连锁,那么我们去掉其中一个突变的所有信息,并不会影响总的信息含量。主成分分析也即是基于这样一种思想开展的,将变量之间根据相关性进行分解、合并和降维,类似于从n维空间到r维空间的投影。
PCA优缺点
优点:
1、以方差衡量信息的无监督学习,不受样本标签限制。
2、由于协方差矩阵对称,因此k个特征向量之间两两正交,也就是各主成分之间正交,正交就肯定线性不相关,可消除原始数据成分间的相互影响
3. 可减少指标选择的工作量
4.用少数指标代替多数指标,利用PCA降维是最常用的算法
5. 计算方法简单,易于在计算机上实现。
缺点:
1、主成分解释其含义往往具有一定的模糊性,不如原始样本完整
2、贡献率小的主成分往往可能含有对样本差异的重要信息,也就是可能对于区分样本的类别(标签)更有用
3、特征值矩阵的正交向量空间是否唯一有待讨论
4、无监督学习
PCA结果解读
在分析过程中,PCA可以让我们非常直观地看出各个样本之间的相似性。例如在一张PCA散点图中,数个样本的点聚在一起,那么就说明这几个样本之间的相似性非常高;反之,如果几个样本的点非常分散,则说明这几个样本之间的相似性比较低。
例1:
下图:几个组的样本对应的散点在组内呈现相互聚集的情况,说明组内的重复性比较好,样本数据非常相似,而组间则有较好的区分度。有的时候为了说明组内样本的相似程度,还会用一个椭圆将同一组的样本对应的散点全部囊括起来。
例2:
并不是所有的PCA结果都是这种明显的分群结果会比较好,还是要根据分析的目的而定。例如在GWAS分析当中,这种“天女散花”一般的PCA散点图,正说明了样本之间不具备明显的亚群分化,适宜进行后续的GWAS分析。
例3:
下图:通过主成分分析方法(PCA)分析9种食物的蛋白质消耗量(变量)与25个欧洲国家(样本)之间的关系
样本点连线距离长 = 样本之间差异性大
样本点连线距离短 = 样本之间差异性小
由图可得,大部分欧洲国家蛋白摄入习惯是:吃鸡蛋、红肉(猪牛羊等畜肉)、白肉(禽、鱼肉及水产品),喝牛奶。
1、各样本点连线的距离:体现各国家蛋白摄入习惯的相似性。
2、主成分与原变量之间的关系:箭头对应的原始变量在投影到水平和垂直方向上后的值,可以分别体现该变量与PC1和PC2的相关性(正负相关性及其大小)(例如,Eggs对PC1具有较大的贡献,而Nuts则与PC1之间呈较大的负相关性)。
3、样本点和箭头之间的距离:反映样本与原始变量的关系。(对于图中用蓝色粗箭头所指的样本点而言,该国的蛋白质来源主要为Fruits and Vegetables)。
例4:
通过主成分分析方法(PCA)分析9种食物的蛋白质消耗量(变量)与25个欧洲国家(样本)之间的关系(3D图)
例3的三维图(取三个主成分的PA图),读图方法不变,只不过是分析样本点在三维空间的距离,而非平面。
例5:
下图:以肠道菌群种类为变量分析肠道菌群相似性
样本间有一定程度的聚类,CC组和HC组的肠道菌群存在明显分离,重叠部分不大(即两个种群相似度不大)。由图可得,结肠癌后肠道菌群发生了一定程度变化。
例6:
下图:以多种蛋白质的表达量为指标区分样本
30例健康组样本点(Healthy)聚集在右上角区域,与75例HCC组和LC组的样本点明显区分。而且HCC组的7个样本点和LC组的2个样本点位于健康者区域,没有健康者被归类为HCC或LC组。
PCA图绘制,R语言代码
数据来源
从TCGA下载结肠癌+直肠癌的表达数据并进行合并,经过配对处理,挑选出44个病人的44个肿瘤样本与44个正常样本。
原始数据"exp_pair"使用count格式,然后经过 dat <- log(exp_pair + 1) 处理形成下图中的数据"dat"。 前44个是正常样本,后44个是肿瘤样本。
dat <- log(exp_pair + 1) # 因为原始数据呈泊松分布,无法画PCA图,故取其log值。
PCA分析函数
做PCA的函数有很多,常用的函数有:princomp,prcomp及rda。
我们在R中输入的数据类型有两类,分别为R mode和Q mode。一般来说数据每一列为一个变量(variable),每一行为一个数据(observation)。其中R mode的数据列数小于行数,是基于列变量的分析;Q mode数据行数小于列数,是基于行数据的分析。如上述表达数据矩阵中样本数小于基因数,属于R mode型数据。
princomp只能用于R mode,它基于协方差(covariance) 或者相关矩阵(correlation) 提取的特征(eigen)并进行特征值分解。默认用法为 x.princomp = princomp(x, cor = FALSE, scores = TRUE)。
Cor是逻辑值的参数,默认cor = FALSE用协方差矩阵计算。cor = TRUE就会用相关矩阵计算特征值。
Princomp 函数的说明文档中推荐 prcomp 更好: The calculation is done using eigen on the correlation or covariance matrix, as determined by cor.
This is done for compatibility with the S-PLUS result.
A preferred method of calculation is to use svd on x, as is done in prcomp.
默认用法为 x.prcomp = prcomp(x, center = TRUE, scale. = FALSE)。
Center为逻辑值,决定是否要以0为中心化,Scale为逻辑值,决定是否要以单位方差进行标准化。
Prcomp 函数的说明文档中说这种方法在数值上更准确:
This is generally the preferred method fornumerical accuracy。
Princomp和prcomp在算法上略有差异。除了分别为特征值分解和奇异值分解外,两者在之前计算协方差的时候标准化的过程存在差异:princomp 计算时分母为N,而 prcomp 分母为N-1。
冗余分析(redundancy analysis,RDA)是一种回归分析结合主成分分析的排序方法,也是多响应变量(multiresponse)回归分析的拓展。从概念上讲,RDA是响应变量矩阵与解释变量之间多元多重线性回归的拟合值矩阵的PCA分析。rda是vegan包的一个函数,可以用rda这个函数来做PCA。虽然简单,但是功能强大。只输入OTU表时做PCA,如果再加上环境因子就做RDA。函数的说明文档中没有专门提做PCA时的方法。但是做RDA采用的是奇异值分解。
princomp和prcomp都是R自带的stats包中的函数。一般用 prcomp 函数即可。
其次,使用 FactoMineR::PCA 也可以做PCA分析。
1. 使用 stats::prcomp 函数进行PCA分析
# 数据转置,转换成行为样本,列为基因的矩阵。prcomp函数要求分析矩阵的列包含特征向量。 data <- t(dat)
# 使用stats::prcomp函数进行PCA分析 prcomp.pca <- stats::prcomp(data) # 查看PCA的分析结果 summary(prcomp.pca) # 结果中PC的数量为nrow和ncol两数中比较小的那个数 # Standard deviation:标准差,其平方为方差 = 特征值 # Proportion of Variance:方差贡献率 # Cumulative Proportion:方差累计贡献率
# 绘制主成分的碎石图 screeplot(prcomp.pca, npcs = 10, type = "lines") # 折线图 screeplot(prcomp.pca, npcs = 6, type = "barplot") # 条形图
# 查看主成分的结果 head(prcomp.pca$x)
2. 使用 FactoMineR::PCA 函数进行PCA分析
# 数据转置,转换成行为样本,列为基因的矩阵。PCA 函数要求分析矩阵的列包含特征向量。 data <- t(dat)
转置结果同上。
# 使用 FactoMineR::PCA 函数进行PCA分析 PCA.pca <- FactoMineR::PCA(data, graph = FALSE) # 查看PCA的分析结果 summary(PCA.pca) # 结果中DIM的数量为nrow和ncol两数中比较小的那个数-1 # Variance:方差, = 特征值 # % of var.:方差贡献率 # Cumulative % of var.:方差累计贡献率 # Dim:Dimension,维度,类似于主成分PC
除了主成分的结果之外,还提供了对于“个体”和“变量”的汇总,即 Individuals,Variables。
Individuals 中,Dist:distance。测量每个个体在多维空间中从原点到样本点的距离。Dim.X:Dimension。测量每个个体在多维空间中从原点到主成分X的距离投影。
ctr:contribution。以百分比形式表示每个个体对给定主成分的贡献。详见 PCA.pca$ind$contrib,每列总和为100%。
cos2:是每个主成分的平方余弦,计算为(Dim.X/Dist)^2。对于给定的主成分,它越接近1,该主成分在捕捉该个体的所有特征方面就越好。
Variables 中,“Dim.X”/“ctr”/“cos2”的解释类似。
# 绘制主成分的碎石图: # FactoMineR::PCA 函数分析的结果用screeplot函数画不出碎石图,总是报错,具体原因不详。但是可以用 factoextra::fviz_screeplot 函数画出来: factoextra::fviz_screeplot(PCA.pca, addlabels = TRUE)
# 查看主成分的结果 head(PCA.pca[["ind"]][["coord"]]) # 结果只显示前5个Dim
3. 使用 graphics 包绘制PCA图
# 使用 graphics 包绘制PCA图 data <- t(dat) # 转置,转换成行为样本,列为基因的矩阵 # 使用stats::prcomp函数进行PCA分析 prcomp.pca <- prcomp(data) # 使用FactoMineR::PCA函数进行PCA分析 PCA.pca <- FactoMineR::PCA(data, graph = F) # 绘制主成分的碎石图 screeplot(prcomp.pca, npcs = 10, type = "lines") factoextra::fviz_screeplot(PCA.pca, addlabels = TRUE) # 绘制stats::prcomp函数分析的PCA图 plot(prcomp.pca$x, cex = 2, main = "PCA analysis", # cex:点大小。main:标题。prcomp.pca$x是 matrix 形式,plot函数会自动取 matrix 形式数据的前两列,即PC1和PC2。 col = c(rep("blue",44), rep("red",44)), # 点颜色,重复次数根据样本数量改动。这里正常样本 pch = c(rep(16,44), rep(17,44))) # 点形状,见?graphics::points。重复次数根据样本数量改动。 # 绘制FactoMineR::PCA函数分析的PCA图,只修改输入数据即可 plot(PCA.pca[["ind"]][["coord"]], cex = 2, main = "PCA analysis", # cex:点大小。main:标题。prcomp.pca$x是 matrix 形式,plot函数会自动取 matrix 形式数据的前两列,即PC1和PC2。 col = c(rep("blue",44), rep("red",44)), # 点颜色,重复次数根据样本数量改动。这里正常样本 pch = c(rep(16,44), rep(17,44))) # 点形状,见?graphics::points。重复次数根据样本数量改动。 # 添加分隔线 abline(h = 0, v = 0, # h:horizontal line,水平线的y值。v:vertical line,垂直线的x值。 lty = 2, col = "gray") # lty:line type,线类型。2表示虚线,详见graphics::par。col:color,线颜色。 # 给每个点添加标签 text(prcomp.pca$x, labels = rownames(prcomp.pca$x), pos = 4, # pos:position,1、2、3和4 分别表示指定(x,y)坐标下、左、上和右侧的位置。 offset = 0.6, cex = 0.5) # offset:偏移量。cex:character expansion,字符大小。 # 添加图例 legend("topleft", title = "Sample", inset = 0, cex = 1, # inset:图例距离图中心的距离。cex:图例字体大小 legend = c("normal", "tumor"), col = c("blue", "red"), pch = c(16, 17))
prcomp函数分析画的PCA图横纵坐标为PC1、PC2。PCA函数分析画的PCA图横纵坐标为Dim1、Dim2。
但是该基础包的函数画出的图比较粗糙。
4. 使用 ggplot2 包绘制PCA图
# 使用ggplot2包绘制PCA图 data <- t(dat) # 转置,转换成行为样本,列为基因的矩阵 library(ggplot2) color = c("#2874C5", "#f87669", "#e6b707", "#868686", "#92C5DE", "#F4A582", "#66C2A5", "#FC8D62", "#8DA0CB", "#E78AC3", "#A6D854", "#FFD92F", "#E5C494", "#B3B3B3") # 设置颜色备选项 # 使用stats::prcomp函数进行PCA分析 prcomp.pca <- prcomp(data) # 使用FactoMineR::PCA函数进行PCA分析 PCA.pca <- FactoMineR::PCA(data, graph = F) # 绘制主成分的碎石图 screeplot(prcomp.pca, npcs = 6, type = "barplot") factoextra::fviz_screeplot(PCA.pca, addlabels = TRUE) # 绘制PCA图 group_list <- as.factor(c(rep(0,44),rep(1,44))) # 设置分组变量,1为实验组,0为对照组 pdat = data.frame(PCA.pca[["ind"]][["coord"]], Group = group_list) # 准备画图输入数据 p = ggplot(pdat, aes(Dim.1, Dim.2)) + # 画背景、x轴、y轴 geom_point(aes(Dim.1, Dim.2, fill = Group), shape = 21, colour = "black") + # 画点图 scale_color_manual(values = color[1:nlevels(group_list)]) + # 改变点颜色 scale_fill_manual(values = color[1:nlevels(group_list)]) + # 改变点填充颜色 theme_classic() + # 修改外观主题,此处为经典主题,即具有x、y轴线,没有网格线 theme(legend.position = "top") + # 调整图例位置,默认在右侧 labs(fill= "Group") + # 修改图例标题,默认为分组变量的列名 stat_ellipse(aes(color = Group, fill = Group), geom = "polygon", alpha = 0.3, linetype = 2) # 添加置信椭圆线,color=线条颜色,fill=填充颜色,alpha=椭圆透明度,linetype=线类型 return(p)
5. 使用 factoextra 包绘制PCA图
# 使用 factoextra 包绘制PCA图 data <- t(dat) # 转置,转换成行为样本,列为基因的矩阵 library(factoextra) color = c("#2874C5", "#f87669", "#e6b707", "#868686", "#92C5DE", "#F4A582", "#66C2A5", "#FC8D62", "#8DA0CB", "#E78AC3", "#A6D854", "#FFD92F", "#E5C494", "#B3B3B3") # 设置颜色备选项 # 使用stats::prcomp函数进行PCA分析 prcomp.pca <- prcomp(data) # 使用FactoMineR::PCA函数进行PCA分析 PCA.pca <- FactoMineR::PCA(data, graph = F) # 绘制主成分的碎石图 screeplot(prcomp.pca, npcs = 6, type = "barplot") factoextra::fviz_screeplot(PCA.pca, addlabels = TRUE) # 添加百分比标签 # 绘制PCA图 group_list <- as.factor(c(rep(0,44),rep(1,44))) # 设置分组变量,1为实验组,0为对照组 # 样本按group_list分组的点图 factoextra::fviz_pca_ind(PCA.pca, geom.ind = "point", # 样本只展示点图 col.ind = group_list, # 样本颜色按group_list区分 palette = color, # 按组着色或填充的调色板 addEllipses = TRUE, # 添加置信椭圆 legend.title = "Groups") # 图例标题 factoextra::fviz_pca_ind(PCA.pca, label = "none", # 文本指定要标记的元素,ind标记样本名,var标记变量名,all标记样本和变量,none不做标记 habillage = group_list, # 按组对观测值进行着色 addEllipses = TRUE, ellipse.level = 0.95, # 椭圆置信度 palette = c("#999999", "#E69F00", "#56B4E9")) # 样本按cos^2分组的点图 # 对于给定的主成分,cos^2越接近1,该主成分在捕捉该个体的所有特征方面就越好。 factoextra::fviz_pca_ind(PCA.pca, geom = "point", # 样本只展示点图 col.ind ="cos2", # 样本颜色按cos^2区分 gradient.cols = c("white", "#2E9FDF", "#FC4E07")) # 梯度颜色 # 变量图,展示各变量箭头,代表主成分载荷 factoextra::fviz_pca_var(PCA.pca, col.var = "steelblue") # 变量图,以百分比形式展示各变量对给定主成分的贡献 factoextra::fviz_pca_var(PCA.pca, col.var = "contrib", # 按贡献度大小区分颜色 gradient.cols = c("white", "blue", "red"), # 梯度颜色 ggtheme = theme_minimal()) # 没有背景注释的简约主题 # 样本和变量的双标图 factoextra::fviz_pca_biplot(PCA.pca, label = "var", # lable:文本指定要标记的元素,ind标记样本名,var标记变量名,all标记样本和变量,none不做标记 habillage = group_list, # 按组对观测值进行着色 addEllipses = TRUE, ellipse.level = 0.95, # 椭圆置信度 ggtheme = theme_minimal()) # 没有背景注释的简约主题
样本按group_list分组的点图:
对于给定的主成分,cos^2越接近1,该主成分在捕捉该个体的所有特征方面就越好。
变量图,展示各变量箭头,代表主成分载荷:
变量图,以百分比形式展示各变量对给定主成分的贡献:
样本和变量的双标图:
6. 使用 scatterplot3d 包绘制PCA图
# 使用 scatterplot3d 包绘制PCA图 data <- t(dat) # 转置,转换成行为样本,列为基因的矩阵 library(scatterplot3d) color = c("#2874C5", "#f87669", "#e6b707", "#868686", "#92C5DE", "#F4A582", "#66C2A5", "#FC8D62", "#8DA0CB", "#E78AC3", "#A6D854", "#FFD92F", "#E5C494", "#B3B3B3") # 使用stats::prcomp函数进行PCA分析 prcomp.pca <- prcomp(data) # 使用FactoMineR::PCA函数进行PCA分析 PCA.pca <- FactoMineR::PCA(data, graph = F) # 绘制主成分的碎石图 screeplot(prcomp.pca, npcs = 6, type = "barplot") factoextra::fviz_screeplot(PCA.pca, addlabels = TRUE) # 添加百分比标签 group_list <- as.factor(c(rep(0,44),rep(1,44))) # 设置分组变量,1为实验组,0为对照组 colors = color[as.numeric(group_list)] # 为每个样本点设置颜色 pdat = data.frame(PCA.pca[["ind"]][["coord"]], Group = group_list) # 准备画图输入数据 # 画PCA函数的点图 scatterplot3d::scatterplot3d(pdat[, 1:3], main = "PCA", # 取前3位主成分,main:标题 pch = 21, cex.symbols = 1.5, color = "black", bg = colors) # pch:点形状。color:点外圈颜色。bg:点填充颜色 # 给PCA函数的点图添加图例 graphics::legend("topright", # 图例位置,有"bottomright", "bottom", "bottomleft", "left", "topleft", "top", "topright", "right" and "center" legend = c("normal", "tumor"), # 图例文本信息来源 pch = 21, col = "black", pt.bg = color[1:nlevels(group_list)], # pch:点形状。col:点外圈颜色。pt.bg:点填充颜色 xpd = TRUE, # FALSE为不允许在作图区域外作图,改为TRUE即可 horiz = FALSE, # 图例放置方向,TRUE为水平放置图例,FALSE为垂直放置 inset = -0.09, # 图例距离图中心的距离 bty = "o") # 图例框是否画出,"o"为画出,"n"为不画出 # 画prcomp函数的点图 scatterplot3d::scatterplot3d(prcomp.pca$x[,1:3], main = "3D PCA plot", # 取前3位主成分。main:标题 pch = c(rep(16,44), rep(17,44)), # 点形状 color = c(rep("blue",44), rep("red",44)), # 点颜色 angle = 45, # angle:观察角度 cex.symbols = 1.5, # cex.symbols:点放大倍率 mar = c(5,4,4,5)) # c(bottom, left, top, right)的数值向量,设置立方体底面、左面、上面、右面到图边缘的距离。 # 给prcomp函数的点图添加图例 graphics::legend("topright", # 图例位置,有"bottomright", "bottom", "bottomleft", "left", "topleft", "top", "topright", "right" and "center" title = "Sample", # 图例标题 xpd = TRUE, # FALSE为不允许在作图区域外作图,改为TRUE即可 inset= -0.11, # 图例距离图中心的距离 legend = c("normal","tumor"), # 图例文本信息来源 col = c("blue","red"), # 图例点颜色 pch = c(16,17)) # 点形状