第1个回答 2022-06-30
title: DALS021-MDS与PCA
date: 2019-08-21 12:0:00
type: "tags"
tags:
这一部分是《Data Analysis for the life sciences》的第8章统计模型的第2小节,这一部分的主要内容涉及MDS和PCA,相应的Rmarkdown文档可以参考作者的 Github 。
MDS的全称为multi-dimensional scaling,即多维数据缩放。在这 一部分中,我们会使用基因表达的数据来作为案例讲解一下。为了简化说明,我们仅考虑3个组织:
结果如下所示:
现在我们要研究一下这个数据集,我们想知道,储存在 mat 列中的基因表达谱的数据在不同的组织间的相似性如何。由于数据很大,无法直接画出相应的多维点图。我们通常只能绘制出二维图形,如果我们要绘制出每两个样本之间的基因表达情况不现实。而MDS图形就是为了解决这个问题而提出来的。
前面我们已经知道了SVD和矩阵代数,那么我们理解MDS就相对清楚了。为了说明MDS,我们先来看一下SVD分解,如下所示:
我们假设 的前两列的平方和剩余列的平方和。因此它们可以写为 其中 是 是第i列(原文是i-th entry)。当出现这种情况时,我们就会得到如下公式:
这就表明,第 列近似等于:
如果我们们定义下面的二维向量:
那么:
上面的这个推导告诉我们,在样本 和样本 之最的距离近拟等于下面二维数据点的距离:
因为 是一个二维向量,因此我们可以通过绘制 和 来发展示这两个样本的距离。现在我们绘制出它们的距离:
从图片上我们可以看出,数据点按照相应的组织区分开来了。上面的这种分开的精确近似取决于前两个主成分解释变异的程度。像上面那样所示,我们可以绘制出每个主成分可以解释的变异程度:
虽然前两个主成分解释了超过50%的变异,不过前面的图形还是没有展示出大量的信息。但是这种图已经足够用于进行可视化大量的数据了。此外,我们还可以注意到,我们能够绘制其它的主成分来研究这些数据点,例如我们绘制第3个和第4个主成分:
从上面图形中我们可以看到,第4个主成分能够将肾脏组织的样本强烈分开。在后面的部分中,我们会讲到批次效应(batch effects)会解释这种情况。
我们在上面使用了 svd() 函数来进行计算,不过R中有一个专门的函数用于计算MDS,生成MDS图。这个函数就是 cmdscale() 函数,这个函数将距离对象作为参数,然后使用主成分分析来对这些距离进行近似计算。这个函数比使用 svd() 函数更高效(因为不可能实现完全的 svd() 函数计算,那样比较花时间)。此函数默认返回二维的数据,不过我们通过设定参数 k (默认情况下, k=2 )可以改变结果中的维度:
再看另外一个:
SVD并非是唯一的,只要我们用 -1 乘以 的样本列,我们就能使用 -1 乘以 的任意列,通过下面的转换我们就能看出来(这一段不懂):
在所有的计算中,当我们计算SVD时,都会扣除行(row)的均值。如果我们要试图计算两列之间的近似距离,那么在 和 之间的距离就与 和 之间的距离相同,因为当我们过计算时,中间的 就会被消去:
因为扣除行均值可以降低总的变异,它可以使得SVD的结果近更为逼近。
P357
PCA的相关资料可以参考作者的 Github 。
前面我们已经提到了PCA,这里继续深入一步,讲一下PCA背后的数学原理。
我们先使用模拟数据的案例展示一个旋转,这个旋转与PCA有着很大的有关系:
这里我们专门来解释一下什么是什么成分(principla components)。
我们使用 这个 矩阵来表示我们的数据。这个类似于我们检测了两组基因的信息,每列表示1个样本。现在我们的任何就是,找到一个 向量 ,使其满足 ,它能使 最大。这个过程可以被视为每个样本,或 向子空间 的投影。因此,我们需要将坐标系进行置换,使新的坐标系能够显示出最大变异。
我先试一下 。这个投影公仅能够给出双胞胎1的身高(橘黄色)。图片标题中显示的是平方和。
我们能否找到一个方向,使得坐标系旋转后,能够表示更高的变异?例如
这个怎么样?它不满足 ,因此我们可以使用另外一个向量,即
这个图形与双胞胎的差异有关,我们知道这个差异很少的。通常平方和我们可以确实这一点,最后我们试一下这个向量:
这个图形与重新缩放(re-scaled)后的平均高度有关,它有着最大的平方和。这是一个数学计算程序,它能够计算出一个 ,能够使平方和最大,SVD就是这样的一个程序。
正交向量能够使平方和最大:
指的就是第1PC。e用于获得PC的加权(weights) 指的就是因子载荷(loadings)。使用旋转这种操作,它指的就是第1PC的旋转方向。
为了获得第2PC,我们可以重复上述操作,但是残差如下:
第2PC的向量含有以下性质:
它能使 最大,
当 是 时,我们可以重复地找到第3,第4,第5,等主成分。
我们已经介绍了如何使用SVD来计算PC。介理,R中有一个专门的函数可以用于找到主成分,即 prcomp() ,在这个案例中,数据默认中心化的,这个函数的使用如下所示:
计算出的结果与SVD相同,直到符号翻转(produces the same results as the SVD up to arbitrary sign flips,实在没理解这句话什么意思)
因子载荷可以通过下面方式计算:
计算结果如下所示:
它就相当于 (up to a sign flip?这个不懂) :
计算结果如下所示:
解释的方差等价于:
计算结果如下所示:
现在我们将 Y 转置一下,因为 prcomp() 函数与我们平时所用的高通量数据储存有点不太一样,平时我们的数据是列为样本,行为特征值,而 prcomp() 函数则是正好相反,它处理的数据列是特征值,行是样本名。