美文网首页
PH525x series - Multi-Dimensiona

PH525x series - Multi-Dimensiona

作者: 3between7 | 来源:发表于2019-12-16 13:50 被阅读0次

多维尺度分析(MDS)

在奇异值分解那一章节我们已将学过:

\mathbf{Y} = \mathbf{U}\mathbf{D}\mathbf{V}^T

在此基础上,我们若假设\mathbf{U}^T\mathbf{Y} = \mathbf{D}V^\mathbf{T}的前两列的平方和远远大于其余列的平方和,换句话说若d_i\mathbf{D}对角线上的第i个奇异值,那么:d_1 + d_2 ≫ d_3 + \cdots + d_n,那么我们便有:
\mathbf{Y}\approx [\mathbf{U}_1 \mathbf{U}_2] \begin{pmatrix} d_{1}&0\\ 0&d_{2}\\ \end{pmatrix} [\mathbf{V}_1 \mathbf{V}_2]^\top
Y_i的第i列便约等于:

\mathbf{Y}_i \approx [\mathbf{U}_1 \mathbf{U}_2] \begin{pmatrix} d_{1}&0\\ 0&d_{2}\\ \end{pmatrix} \begin{pmatrix} v_{i,1}\\ v_{i,2}\\ \end{pmatrix} = [\mathbf{U}_1 \mathbf{U}_2] \begin{pmatrix} d_{1} v_{i,1}\\ d_{2} v_{i,2} \end{pmatrix}

如果我们设Z_i等于:

\mathbf{Z}_i=\begin{pmatrix} d_{1} v_{i,1}\\ d_{2} v_{i,2} \end{pmatrix}那么:
\begin{align*} (\mathbf{Y}_i - \mathbf{Y}_j)^\top(\mathbf{Y}_i - \mathbf{Y}_j) &\approx \left\{ [\mathbf{U}_1 \mathbf{U}_2] (\mathbf{Z}_i-\mathbf{Z}_j) \right\}^\top \left\{[\mathbf{U}_1 \mathbf{U}_2] (\mathbf{Z}_i-\mathbf{Z}_j)\right\}\\ &= (\mathbf{Z}_i-\mathbf{Z}_j)^\top [\mathbf{U}_1 \mathbf{U}_2]^\top [\mathbf{U}_1 \mathbf{U}_2] (\mathbf{Z}_i-\mathbf{Z}_j) \\ &=(\mathbf{Z}_i-\mathbf{Z}_j)^\top(\mathbf{Z}_i-\mathbf{Z}_j)\\ &=(Z_{i,1}-Z_{j,1})^2 + (Z_{i,2}-Z_{j,2})^2 \end{align*}

也就是说,在这种情况下,样本i和样本j之间的距离就可以用两个二维的点进行约算了:
(\mathbf{Y}_i - \mathbf{Y}_j)^\top(\mathbf{Y}_i - \mathbf{Y}_j) \approx (Z_{i,1}-Z_{j,1})^2 + (Z_{i,2}-Z_{j,2})^2

举例说明:

library(rafalib)
library(tissuesGeneExpression)
data(tissuesGeneExpression)
colind <- tissue %in% c("kidney","colon","liver")
mat <- e[,colind]
group <- factor(tissue[colind])

s <- svd(mat-rowMeans(mat))
PC1 <- s$d[1]*s$v[,1]
PC2 <- s$d[2]*s$v[,2]
mypar(1,1)
plot(PC1,PC2,pch=21,bg=as.numeric(group))
legend("bottomright",levels(group),col=seq(along=levels(group)),pch=15,cex=1.5)
201912161109.png

到这一趴的时候我头一开始特别没理解为啥在给这些点涂色的时候直接as.numeric(group)就行,后来想明白了,也就是奇异值分解得到的右奇异矩阵中的每一列(也就是每一个特征向量)中的元素顺序是与原始矩阵的样本顺序一致的。

从上图中我们可以看到单是用前两个维度就可以将这些样本很好的分开,再看下方差解释度的趋势图:

plot(s$d^2/sum(s$d^2))
201912161119.png

可以看到,其实前两个维度的方差解释度加起来也就50%多一些,也就说虽然有很多信息不能在图中展示出来,但是单是用前两个维度画图却仍能够很好的看出点与点之间的分布情况。

我们再看下用第3和第4个主成分做出的图:

PC3 <- s$d[3] * s$v[,3]
PC4 <- s$d[4] * s$v[,4]
mypar(1,1)
plot(PC3,PC4,pch=21,bg=as.numeric(group))
legend("bottomright",levels(group),col=seq(along=levels(group)),pch = 15,cex=1.5)
201912161314.png

在上图中,kidney组织的样本很明显的被分开了,这可能与后面要学的batch effects有关,更加详细的以后再说,关于这点作者在这里也就提了一句。

cmdscale

我们除了可以使用svd()函数画多维尺度图之外,还可以使用cmdscale()函数,它以一个距离矩阵作为输入参数,然后利用PCA为这一矩阵提供一个最近似的k维矩阵(k默认为2),与svd()相比,运算时间更短。

举例说明:

d <- dist(t(mat))
mds <- cmdscale(d)
mypar()
plot(mds[,1],mds[,2],bg=as.numeric(group),pch=21,
     xlab="First dimension",ylab="Second dimension")
legend("bottomleft",levels(group),col=seq(along=levels(group)),pch=15)
201912161330.png

符号问题

观察下图,可以看出这两种方法得到的第一第二主成分的正负符号并不完全一致:

mypar(1,2)
for (i in 1:2){
  plot(mds[,i],s$d[i]*s$v[,i],main=paste("PC",i))
  b <- ifelse(cor(mds[,i],s$v[,i])>0,1,-1)
  abline(0,b)
}
201912161332.png

这是因为奇异值分解的结果并不唯一,关键就在这正负符号上:

\mathbf{-1UD(-1)V}^\top = \mathbf{UDV}^\top

减去平均数的问题

如果你仔细看过文章内容,就会发现在我们之前奇异值分解的计算当中,都会先减去行均值,也就是用原始矩阵减去行均值得到的新矩阵计算距离与用原始矩阵去计算距离是一回事:

\left\{ ( \mathbf{Y}_i- \mathbf{\mu} ) - ( \mathbf{Y}_j - \mathbf{\mu} ) \right\}^\top \left\{ (\mathbf{Y}_i- \mathbf{\mu}) - (\mathbf{Y}_j - \mathbf{\mu} ) \right\} = \left\{ \mathbf{Y}_i- \mathbf{Y}_j \right\}^\top \left\{ \mathbf{Y}_i - \mathbf{Y}_j \right\}

所以,为了减少计算量,我们通常都会这么去处理。

阅读原文请戳

相关文章

  • PH525x series - Multi-Dimensiona

    多维尺度分析(MDS) 在奇异值分解那一章节我们已将学过: 在此基础上,我们若假设的前两列的平方和远远大于其余列的...

  • PH525x series - Exercises - Line

    本篇文章是PH525x series课程中Linear models and randomness的练习章节,下面...

  • 线性回归模型

    在学习PH525x series - Chapter 5 - Linear Models时,觉得有些地方理解起来有...

  • PH525x series - Hierarchical Mod

    在上一篇文章PH525x series - Bayesian Statistics中是将层次模型应用到了棒球运动当...

  • PH525x series - Collinearity

    共线性 当自变量之间存在共线性时,线性回归得到的最小二乘估计的值并不唯一。共线性简单点说就是,设计矩阵中的某几列存...

  • PH525x series - Introduction to

    本章会对线性模型做一个大致的介绍,还是举例说明吧: 例1:自由落体问题 想象自己是16世纪的伽利略,正在研究自由落...

  • PH525x series - Projections

    前面的章节学的是降维、奇异值分解以及主成分分析的大致内容,本篇文章则开始更加详细的介绍这背后的数学原理,首先要学的...

  • PH525x series - Running PCA and

    在PCA相关的章节最后,系列教程的作者又专门写了一章“在R中运行PCA和SVD”,使用的还是tissuesGene...

  • PH525x series - Statistical Mode

    正连续值的分布 在生物学中有很多数据的分布特征是“strictly positive and heavy righ...

  • PH525x series - Principal Compon

    这一章,作者就是在数学原理方面又细讲了下主成分分析(PCA) 例子:双胞胎身高 作者首先使用双胞胎身高的例子来说明...

网友评论

      本文标题:PH525x series - Multi-Dimensiona

      本文链接:https://www.haomeiwen.com/subject/womrnctx.html