03高通量测序-PCA中的主要概念

如题所述

第1个回答  2022-06-29

我们使用SVD(singular value decomposition,中文译名“奇异值分解”)的方法来计算PCA。我们先看一个简单的案例,在这个案例中,我们检测了6只不同小鼠的2个基因。其实我们可以把它再抽象化一下,把小鼠看成样本,基因看成2个变量。如果我们只检测1个基因的话(Gene 1),那么我们根据基因1表达的情况,把小鼠的绘制到数轴上,小鼠1,小鼠2,小鼠3的Gene 1表达水平比较高,而小鼠4,小鼠5,小鼠6的Gene 1水平则较低。虽然这个图形比较简单,但是,我们从中还是能得到一些信息的,例如小鼠1,小鼠2和小鼠3比较接近,小鼠4,小鼠5和小鼠6比较接近。

如果我们检测了2个基因,那么我们可以绘制一个二维坐标系,横轴是Gene 1,纵轴是Gene 2,那么小鼠1,小鼠2和小鼠3会聚在一起,小鼠4,小鼠5和小鼠6会聚在一起。

如果我们检测了3个基因,那么我们可以绘制三维的坐标系,在上图的这个3维坐系中,圆点越大,表示离你越近。如下所示:

再进一步,如果我们检测了4个基因,此时我们很难绘制出四维的坐标系,那么我们就需要进行PCA分析了,PCA可以把超过4个的基因降维成二维的坐标系,在这个PCA的二维坐标系中,我们可以发现,小鼠4、小鼠5和小鼠6是一类,小鼠1,小鼠2和小鼠3是一类,它们的各自的基因表达模式也类似,PCA在对数据进行聚类(clustering)时有很大的价值,例如,经过PCA分析,在它的二维坐标轴上我们可以发现,Gene 3在x轴上对样本的区分有贡献最大

我们还回最初2个基因的案例上来:

此时,我们可能有疑问,为什么这条直线是最匹配数据的,它的计算原理是什么,那么接着看。我们先回到最初的直线,为了准确地找出最佳匹配所有数据的直线,PCA会将所有数据点都映射到这条直线上来,此时, 可以计算这些数据点到投射到这条直线上的距离,并且使这些距离最小,除了可以计算数据点到直线的距离最小外,还要计算所有数据点投射到这条直线上的点(图中绿叉位置所在点)到原点的距离,使这个距离最大。

通过勾股定理,我们可以知道,因为a不变,当数据点到直线的距离最短时 (b) ,投影点到原点的距离最大 (c)

我们计算投影点到原点的距离,我们把它我们把它命名为d1,d2......计算剩余的投影点到原点的距离。然后把这些值的平方加起来称为SS(distances)。我们旋转直线,直到SS的值最大,此时数据点到直线的距离最短。这条直线就叫第一主成分(Principal Component 1,简称PC1)。

对于PC1,斜率为0.25。也就是说基因1增加4个单位,基因2增加1个单位。计算出红色箭头的长度为4.12。

当你用SVD(singular value decomposition,奇异值分解)进行PCA时,红色箭头的长度=1,我们所要做的就是把这个三角形缩小到红箭头是1个单位时,只需每边除以4.12。

三个边长分别变成了1,0.242,0.97,但,Gene1/Gene2仍然等于4。此时我们可以说PC1由0.97的Gene1和0.242的Gene2构成。

我们回顾一下计算过程:

这个单位向量由基因1的0.97和基因2的0.242部分组成, 称为PC1的“奇异向量”(Singular Vector)或“特征向量”(Eigenvector) ,每个基因的比例则被称为 载荷得分(Loading Scores) 。原始数据的投影点到原点的距离的平方SS被称为PC1的 特征值(Eigenvalue) 。PC1的特征值的平方根叫做PC1的 奇异值(Singular value)

因为是一个二维图, PC2 只是一条垂直于PC1的穿过原点的直线,没有任何进一步的优化要做。由于PC2与PC1垂直,所以斜率为-4,也就是说PC2由-1份Gene1和4份Gene2组成。如果我们对所有东西进行缩放,得到一个单位向量,PC2由-0.242个Gene 1和0.97个Gene 2构成, 称为PC2的“奇异向量”(Singular Vector)或“特征向量”(Eigenvector) 。每个基因的比例则被称为 载荷得分(Loading Scores) ,告诉我们,就基因值如何投射到PC2上而言,Gene2的重要性是Gene1的4倍

最后,PC2的 特征值 是投影点到原点的距离的平方和。

此时,PC1和PC2的计算结束,绘制最终的PCA图,如下所示:

然后旋转这个坐标,让PC1水平,PC2垂直,如下所示:

在这个新的坐标系中,图中黑色的叉就表示原始的样本6(Sample 6),如下所示:

而Sample 6位于这个点上:

同理,Sample 2在这里:

我们可以将特征值转化为PC1到原点的变异,通过除以样本大小减1:n-1。这个例子,假设PC1的变异为15,PC2的变异为3。这意味着PCs的变化是15 +3= 18。PC1占了PCs变异的15 / 18= 0.83 = 83%。PC2占了3/18= 0.17= 17%的PCs变异。

碎石图(scree plot) 是用图形表示每个PC所占的变异百分比。

PCA有3个变量(在这种情况下,3个基因)几乎等同于2个变量

然后找到经过原点的最佳拟合直线,和之前一样,最佳拟合线是PC1。PC1现在有三种成分:0.62的Gene1、0.15的Gene2和 0.77的Gene3,Gene3 是最主要的组成部分。然后求出PC2,它经过原点并垂直于PC1。PC2现在有三种成分:0.77的Gene1、0.62的Gene2和 0.15的Gene3,Gene1 是最主要的组成部分。然后,我们找到了PC3,这条最合适的直线,它通过原点并垂直于PC1和PC2。如果我们有更多的基因,我们就会通过添加垂线和旋转它们来不断寻找越来越多的主成分, 理论上,每个基因(或变量)都有一个。但在实际操作中,PC的数量不是变量的数量或者样本的数量,取其中较小的一个。 一旦你找出了所有的主成分,你可以使用特征值(即SS(距离))来确定每个PC的变化比例。在这个例子中,PC1=79%,PC2=15%,PC3=6% ,PC1和PC2占了变异的绝大比例。 这就表明了,在二维图中,我们基本上只使用PC1和PC2就能解释三维图中的数据,因为二维图中的PC1和PC2占据了整体的变异的94%,

为了将3-D图像转换成2-D的PCA图像,我们去掉了所有除了数据和PC1、PC2。将数据投影到PC1和PC2上,然后旋转坐标轴,这是我们新的PCA图中的样本4。

最后,我们使用PC1和PC2将数据绘制成二维图形。

如果我们测量每只老鼠的4个基因,我们不可能画出一个四维图数据,但这并不妨碍我们进行PCA计算,并查看碎石图。在这种情况下,我们可以计算主成分,发现PC1和PC2占变异的90%,所以我们可以使用它们来绘制二维PCA图。

注意:如果碎石图中PC3和PC4占据了大量的变化,那么仅仅使用前2个PCs并不能创建一个非常准确的数据表示。然而,即使像这样一个PCA图也可以用来识别数据分类。