查看: 51132|回复: 106

【主题帖】Omicshare tools——二维PCA分析使用教程

  [复制链接]

管理员

Rank: 15Rank: 15Rank: 15Rank: 15Rank: 15

主题
49
注册时间
2015.12.5
在线时间
533 小时

活跃会员论坛元老


发表于 2016.4.27 02:30:50 | 显示全部楼层 |阅读模式
PCA分析是一种常见的样本聚类的方法,常用于基因表达、16s 多样性分析、重测序等基于各类变量信息(各个基因表达量、各个菌落丰度、SNP的基因型值)的样本聚类,以下就将针对OS-tools中的”主成分分析(二维PCA图)“进行介绍。
一、数据准备——PCA分析的输入文件
必填项:“上传矩阵表格的文件”
  这个选项你需要上传典型的表达量表格文件,必须是TXT格式的文本文件。在这个文件中,表头(第一行)是样本名称;第一列是基因名称,表中其他元素是基因表达量值。为了保证输入文件没有错误,建议阅读以下的帖子:
OS-tools 常见使用错误(不定期更新)
(出处: OmicShare Forum)

选填项目ID文件
    这个文件可以改变PCA分析涉及的样本和分组信息,有两列:
    第一列是样本的ID,就是填写你需要用于PCA分析的样本。如果是所有样本都纳入分析,当然就写所有样本的ID了。但如果进行所有样本的PCA分析后,你发些有离群样本。由此你打算剔除离群样本重新评估其他样本的聚类效果。那么你就可以重新选择其他样本进行PCA分析。
    第二列是分组信息。如同范例文件演示的,一共有6个样本,分为两组,分别是C组(control组)和M组(mutant组)。那么,你就可以在第二列写入分组信息。那么进行PCA分析的时候,软件将可以根据分组信息,在PCA分析的结果图中将不同组的样本打上不同的颜色,结果更加简单易懂,结果如图1。
    当然,如果你没有分组的信息的化,第二列依然可以写上样本ID,那么结果就是1个样本被定义为1组。
   如果你不提供这个文件,软件将使用所有样本进行分析,1个样本就代表1组。
选填项目ID文件
   输入文件也是1个文本文件,列出用于PCA分析的基因列表。如果你不提供这个文件,我们将用所有基因用于PCA分析。这个参数的用途可能有两个情境:
    a)  你发些样本中某些基因表达异常,甚至是外源污染,其影响了聚类效果。那么,你可以提供1个不含异常基因的列表,进行PCA分析;
    b)  你不想使用所有基因,而是某些目标基因(差异表达基因,某些通路或某些功能的基因)进行PCA聚类,看这些基因对样本进行分类的效果,那么就可以选择这个参数进行分析。
二、PCA结果的解读
PCA分析结果的打包文件中,对用户有用是PCA结果的聚类散点图。比如,图1是拟南芥两组样本的PCA分析结果。图中使用主成分1(PC1)和主成分2作为X轴和Y轴绘制散点图,每个点代表1个样本。在这样的PCA图中,如果两个样本距离越远,则说明两个样本基因的表达模式差异越大。反之,则说明相应样本整体表达模式越接近。所以,PCA分析常用于评估样本重复性的好坏。理想情况下,生物学重复的样本应该聚类在一起,而处理组间应该可以清晰区分开。如图1 中,在PC1的维度上,C组和M组可以清晰地区分开,所以总体而言两组实验组的区分度还是很高的。但C2样本看起来离开C1、C3很远,是否代表C2样本的重复性不好呢?
    回答这个问题,你需要注意在在X轴和Y轴的PC1和PC2后的括号中有一个百分比。这个百分比代表这个PC1或PC2这个两个变量可以解释的总体方差的百分比。PC1占到总体方差变异的99.7%,所以评估这6个样本的分类关系,只要看PC1就足够了。PC2占的总体方差的信息量极少(0.2%),可以忽略不计。所以,C2虽然离群,在因为与C1、C3的差异(距离)主要都在PC2这个维度上,这个影响可以忽略不计。而C组和M组的差异,主要在PC1这个维度上,说明处理组间的差异极大,大大超过组内的差异,总体而言6个样本分类效果还是不错的,说明实验结果还是很成功的。
    看到这里,你基本可以对PCA分析的结果图进行解读了,那你也可以到此为止了。
    但你可能还是有一些疑问,PC1是什么东西?是1个基因吗?占到总体方差的百分比是什么意思?这个就涉及到PCA原理的解释。
图1 PCA结果
三、PCA分析的原理详解
    PCA分析是一种线性代数的处理方法,我们可以从一个更为直观的空间立体几何的角度来解释这个问题。
    PCA分析应用的情境就是:当一组样本有大量指标(例如:10000个基因的表达量信息)的情况下,我们如何实现简单直观地对样本进行区分?
    什么是指标?
    如果有一屋子人,要评估哪些人轻,哪些人重,很简单:称量一下每个人的体重就ok了。因为这里只有一个指标,或者说只有1个变量——体重,所以问题非常简单明了。
    如果我们只有两个变量(例如,身高 + 体重)的时候,则依然可以使用二维(X、y轴)散点图将样本区分开。当然,如果增加到3个变量(例如,身高、体重、体脂含量),使用三维直角坐标系(X、Y、X)还是勉强可以胜任的。但如果有10000个变量呢?我们将如何对样本进行分类呢?这就涉及到主要成分的提取以及后续的数据降低维度的过程。
    对于一组样本如果有两个指标(假设为两个基因X1和X2的表达量),将其绘制在二维坐标系中,则如图2(左)。你可以注意到了,如果将样本映射到几何空间结构中后,每个变量(基因)对应着1个维度的坐标轴(X轴和Y轴)。2个基因就组成了1个2维的空间。我们可以注意到,在这组数据的离散程度(总体方差),可以拆分为在X1维度和X2维度上的离散程度。在原来的直角坐标系中,X1和X2方向上这些样本分布范围都是相似的。或者说,X1和X2分别可以解释总体方差的各50%。所以,如果我仅仅保留X1的信息而将X2剔除,虽然变量变少、数据变简单了,但同时也意味着我损失了X2维度的变异信息,即损失了整体变异50%的信息。
    是否有一种,既能降低数据维度又可以尽量减少信息损失的方法呢?比较简单的方法是,在不改变这些点的位置的前提下(样本间的关系没有改变),对坐标轴进行旋转,得到新的坐标轴体系(F1和F2)。在新的坐标轴体系中,变量的数量依然相同(2个),但这些样本在F1和F2方向上的离散度是不同的。F1对总体方差的贡献度更大(图中约为70%),F1含有更多的信息。对应的,F2含有的信息则更少(可能仅仅占总体反差的30%)。那么,现在如果我们将F2丢弃,仅仅保留F1,效果图如图2(右)。尽管只保留了F1的信息,但样本大部分离散关系的信息(约70%)被保留了下来,F1的方向上依然可以较好地对样本进行区分。这就实现了降低简化数据的维度(2维到1维)而信息损失较少(仅30%)的目的。

图2 PCA分析中的主成分提取和降维过程
备注:与《三体》中的降维攻击无关

    大家注意这个过程的特点:
    (1)样本位置并没有变化,变化的是坐标轴;
    (2)坐标轴(变量)的数量没有变化,但不同坐标轴(变量)包含的信息量被重新分配(贫富差距拉大);
    (3) 由于不同变量的“贫富差距”拉大,所以我们可以放弃那些包含信息量较少的变量,从而在损失较小的情况下,实现降低数据的维度(简化数据)。
    (4)这些新的变量(F1和F2),就是对应地被称为主成分1(PC1)和主成分2(PC2。PC1和PC2并非真实的变量(或者说真实的基因表达量),而是所有基因表达量信息共同加权得到的变量(坐标旋转后的夹角其实就是加权系数)。在我们Omicshare tools分析结果的压缩包中的“temp.dat.PC_2Dcomp.xls”文件,就是所有基因对应PC1和PC2的加权系数。

    以上就是直观上的PCA分析过程,当然实际计算过程还是线性代数计算过程,感兴趣的读者可以自己百度相关资料。
    我们可以再略微提高以上情况的复杂程度。例如,我们现在有10个基因,每个基因都占了总体方差的一部分(最高16%,最低7%),如图3(左)。那么,如果我们放弃几个基因的信息,都意味着大量有效信息的丢失。如果我们进行PCA分析,对总体方差的分布进行最大程度的优化重分配,结果如图3(右)。处理后,PC1解释了总体方差的60%,PC2解释了总体方差的32%。PC1+PC2 已经可以解释总体方差的92%。即使我放弃PC3~PC10,将数据的维度从10维降低到2维,也仅仅损失信息量的8%而已。这个损失量还是可以忍受的。既然是二维数据了,问题就简单了,绘制1份二维散点图就可以直观的呈现样本间的聚类关系了(图1的效果)
图3 10个基因的PCA分析过程
    在实际的RNA-seq样本中,变量的数量还不仅仅是10个,而是数以万计的基因。但降维的原理和过程依然是相同的。当然,图1的数据是属于比较极端的情况。因为实际6个样本是拟南芥的野生型组和突变体组比较,处理组内高度相似,组间差异较大,加之进行PCA分析的基因是我们挑选过的.所以PCA分析后,PC1占了总体方差的99.7%,PC2~PCn几乎都可以忽略不计了。在实际的RNA-seq实验数据中,一般PC1+PC2也能占到总体方差的50%以上,所以都可以较好的实现降维(10,000 → 2)的同时,信息损失最小的目的。
四、应用实例
    PCA分析,可以应用在任何 多样本多变量的数据聚类中。基因表达量数据、16s 多样性分析、蛋白定量数据等都可以用OS-tools的PCA分析工具进行分析。当然,其他类型的数据也可以分析,例如生态学调查中多个地点、多个环境指标的PCA分析。重测序的SNP数据,理论上也可以进行PCA分析,但需要先将基因型转化为数字,所以目前OS-tools的PCA工具还不适用于重测序数据的PCA分析。
    PCA分析常用于解释以下几类问题:
    (1)通过PCA散点图中的样本的分布,初步评估实验重复性。在理想条件下,在PCA图中,生物学重复的样本应该聚类在一起;
    (2)观察样本中是否存在离群样本;
    (3)通过PCA样本的分布,判断实验处理效应的大小或实验处理何时起作用。

    PCA的结果,并非1个纯数学问题,应该结合具体的生物学背景对图形进行解读。 以下,我们举两个例子介绍PCA的解读应用:

图4  不同生长时间点(40天和60天)的某植物对不同供水量的应答。
    在图4的PCA图中,我们可以解读出两个信息:
    (1)对于不同的供水量处理,60天的个体分布更加离散。推测是因为60天这个时期属于生成发育的高峰期,个体对供水量的变化更加敏感;
    (2)IL60这个个体属于离群样本,比较异常。可能原因是这个个体遗传背景比较特殊,或实验处理不稳定、或样本污染等。那么在进行深一步的分析前,需要进一步调查这个样本的数据,查明离群的原因,否者可能会影响最终数据的可靠性。


图5 热应激胁迫前后不同时间点的样本RNA-seq数据的PCA结果

(Franssen S U, Gu J, Bergmann N, et al. PNAS,2011, 108(48): 19276-19281.)。

S:south ;N:north;H:heat ;C:control;heat_vave:热应激;Recovery:应激后恢复

  图5来源一篇 PNAS文章,背景是研究一种水生植物 北方种(不耐热)与南方种(耐热)的热应激应答比较。蓝圈中的材料是热应激处理前(以及不处理的对照组)的样本,处于相同的状态,PCA图中聚类在一起。在热应激处理过程中,南北方种在PCA图中的聚类位置发生了移动(红圈中的样本),但依然聚类在一起。从中我们可以推测南北方种的确发了逆境胁迫应答,而且应答反应是相似的(聚类位置一致)。在热应激处理后的恢复期,南北方种的PCA图中的聚类位置却发生了分离。南方种回到了原始的位置(蓝圈),但北方种却依然远离原始的位置(绿圈)。说明,南方种作为耐热品种,在热应激后很快可以恢复到原始状态;但北方种作为不耐热品种,应激造成了难以恢复的细胞内代谢紊乱或损伤,在恢复期依然难以快速恢复到原始状态。这个PCA图形与植物逆境胁迫中的相关规律是一致的:即,高抗和抵抗的植物,在逆境下的状态是类似的,但区别在于高抗性的品系可以在应激结束后迅速恢复到正常生长状态,但低抗性品系却容易发生难以逆转的损伤。

  我们从这两个案例可以看出,PCA结果(包括其他图形)必须结合项目的生物学背景进行解读,才有生物学意义。
五、引用
  由于OS-tools目前还没有发表文章(在我们计划中)。如果大家发表的文章中要引用OS-tools,可以采用以下方法引用:  
  PCA分析的引用方法:
  Principal component analysis(PCA)  was performed using the OmicShare tools,a free online platform for data analysis (www.omicshare.com/tools)

本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有帐号?立即注册

x
新的一天加油!
回复

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
1
注册时间
2016.1.14
在线时间
32 小时

发表于 2016.4.27 08:12:57 | 显示全部楼层
期待以后有更多类型图的解析,这样再外行的人不用过多久也不用再外行了!基迪奥棒棒哒,绝对支持!
回复 支持 1 反对 0

使用道具 举报

版主

Rank: 10Rank: 10Rank: 10

主题
3
注册时间
2016.3.23
在线时间
236 小时

突出贡献优秀版主


发表于 2016.4.27 09:05:39 | 显示全部楼层
真的很好
我要看看,我到底有多么的喜欢你!
回复

使用道具 举报

帝王蝶

Rank: 4

主题
11
注册时间
2016.4.7
在线时间
146 小时

发表于 2016.4.27 15:11:27 | 显示全部楼层
多谢解答
回复

使用道具 举报

草履虫

Rank: 2

主题
0
注册时间
2016.4.22
在线时间
2 小时

发表于 2016.4.27 16:33:47 | 显示全部楼层
写得很好!大写的赞!
回复 支持 1 反对 0

使用道具 举报

版主

Rank: 10Rank: 10Rank: 10

主题
52
注册时间
2016.1.8
在线时间
274 小时

突出贡献优秀版主论坛元老


发表于 2016.4.27 16:59:15 | 显示全部楼层
老师辛苦了,谢谢
回复 支持 反对

使用道具 举报

帝王蝶

Rank: 4

主题
2
注册时间
2016.4.10
在线时间
36 小时

发表于 2016.4.27 17:21:29 | 显示全部楼层
通俗易懂
回复

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
26
注册时间
2016.1.8
在线时间
431 小时

发表于 2016.4.28 14:28:42 | 显示全部楼层
很有用 谢谢
新的一天加油!
回复 支持 反对

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
1
注册时间
2016.4.28
在线时间
2 小时

发表于 2016.4.28 19:51:32 | 显示全部楼层
说的很详细啊
回复 支持 反对

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
10
注册时间
2016.4.7
在线时间
215 小时

活跃会员突出贡献论坛元老


发表于 2016.4.29 10:28:54 | 显示全部楼层
学习了,谢谢!
回复 支持 反对

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
0
注册时间
2016.4.21
在线时间
10 小时

发表于 2016.4.29 13:00:47 | 显示全部楼层
,给力给力,懂了点
回复 支持 反对

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
17
注册时间
2015.12.7
在线时间
145 小时

优秀版主


发表于 2016.4.29 14:35:19 | 显示全部楼层
最后一句亮了
回复 支持 反对

使用道具 举报

功夫熊猫

Rank: 10Rank: 10Rank: 10

主题
1
注册时间
2016.4.25
在线时间
794 小时

灌水之王


发表于 2016.4.29 20:28:34 | 显示全部楼层
学习学习一下
h OK啦我在办公室
回复 支持 反对

使用道具 举报

中华鲟

Rank: 5Rank: 5

主题
9
注册时间
2016.4.12
在线时间
58 小时

发表于 2016.5.1 20:32:02 | 显示全部楼层
解析灰常重要
加油啊
回复 支持 反对

使用道具 举报

帝王蝶

Rank: 4

主题
0
注册时间
2016.3.31
在线时间
129 小时

发表于 2016.5.3 19:30:09 | 显示全部楼层
周老师威武,专业…
回复 支持 反对

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
19
注册时间
2016.4.26
在线时间
250 小时

发表于 2016.5.3 23:47:18 | 显示全部楼层
收藏了!
搞定搞定!!
回复

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
1
注册时间
2016.4.8
在线时间
37 小时

发表于 2016.5.5 15:54:17 | 显示全部楼层

求助,测序pca图是这样,何解?

mirna测序,两个处理,两个重复,pca测出来是这样的



本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有帐号?立即注册

x
回复 支持 反对

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
1
注册时间
2016.4.8
在线时间
37 小时

发表于 2016.5.5 15:56:18 | 显示全部楼层
自己顶一下,捉急,求教大神
回复 支持 反对

使用道具 举报

中华鲟

Rank: 5Rank: 5

主题
14
注册时间
2016.1.4
在线时间
212 小时

发表于 2016.5.6 08:06:35 | 显示全部楼层
红色是一对重复,绿色是一对重复?
回复 支持 反对

使用道具 举报

中华鲟

Rank: 5Rank: 5

主题
14
注册时间
2016.1.4
在线时间
212 小时

发表于 2016.5.6 08:08:12 | 显示全部楼层
你纵坐标占的比重是多少没显示呢?首先横纵坐标的百分比加起来要超过70%呢
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

快速回复 返回顶部 返回列表