获取pheatmap聚类后和标准化后的结果

祥强6csdm0n3vs 2020-09-18

展开全文

前言

NGS系列文章包括NGS基础、转录组分析（Nature重磅综述|关于RNA-seq你想知道的全在这）、ChIP-seq分析（ChIP-seq基本分析流程）、单细胞测序分析 (重磅综述：三万字长文读懂单细胞RNA测序分析的最佳实践教程（原理、代码和评述）)、DNA甲基化分析、重测序分析、GEO数据挖掘（典型医学设计实验GEO数据分析 (step-by-step) - Limma差异分析、火山图、功能富集）等内容。

pheatmap是简单常用的热图绘制包，可以快速、简单、可定制的绘制漂亮热图。具体见R语言学习-热图简化和免费高颜值可定制在线绘图工具 ImageGP。

现在要解决的一个问题是图出来了，想看下转换后用于绘图的表格，也就是获取聚类后的矩阵和聚类标准化后的矩阵。

生成测试数据

mat <- matrix(rnorm(30), nrow=5)

colnames(mat) <- paste('sample', 1:6, sep='_')

rownames(mat) <- paste('gene', 1:5, sep='_')

mat

结果如下

##          sample_1   sample_2   sample_3    sample_4   sample_5   sample_6
## gene_1 -0.3286368  0.3153119 -0.7730821 -0.85242874 -0.5303812  0.5088226
## gene_2 -1.3153020  0.3193550  0.4496518 -1.08782734  1.7620763 -0.9312810
## gene_3  0.6545161 -0.8220414 -1.1916559  0.04775437  0.2814619  1.8720241
## gene_4  1.0810986  0.2298092 -0.3615045  0.70162614  1.8572989  0.7250737
## gene_5 -1.8931573  2.7013864  0.5049798 -0.13541785 -1.7796036 -0.3185864

绘图

library(pheatmap)

# 绘图同时存储绘图结果
(a <- pheatmap(mat, cluster_rows = T, cluster_cols = T))

提取聚类后的原始矩阵

# 查看绘图数据的结构
# 直接查看会很大，这里只展示其前2层
# str: structure
str(a, max.level = 2)

# Rstudio中
# View(a)

结果如下

## List of 4
##  $ tree_row:List of 7
##   ..$ merge      : int [1:4, 1:2] -1 -4 -2 -5 -3 1 2 3
##   ..$ height     : num [1:4] 2.4 3.21 4.38 5.56
##   ..$ order      : int [1:5] 5 2 4 1 3
##   ..$ labels     : chr [1:5] 'gene_1' 'gene_2' 'gene_3' 'gene_4' ...
##   ..$ method     : chr 'complete'
##   ..$ call       : language hclust(d = d, method = method)
##   ..$ dist.method: chr 'euclidean'
##   ..- attr(*, 'class')= chr 'hclust'
##  $ tree_col:List of 7
##   ..$ merge      : int [1:5, 1:2] -1 -6 -2 -5 3 -4 1 -3 2 4
##   ..$ height     : num [1:5] 1.98 2.29 2.55 3.78 5.21
##   ..$ order      : int [1:6] 2 3 5 6 1 4
##   ..$ labels     : chr [1:6] 'sample_1' 'sample_2' 'sample_3' 'sample_4' ...
##   ..$ method     : chr 'complete'
##   ..$ call       : language hclust(d = d, method = method)
##   ..$ dist.method: chr 'euclidean'
##   ..- attr(*, 'class')= chr 'hclust'
##  $ kmeans  : logi NA
##  $ gtable  :List of 6
##   ..$ grobs        :List of 6
##   ..$ layout       :'data.frame':    6 obs. of  7 variables:
##   ..$ widths       :List of 6
##   .. ..- attr(*, 'class')= chr [1:2] 'unit.list' 'unit'
##   ..$ heights      :List of 5
##   .. ..- attr(*, 'class')= chr [1:2] 'unit.list' 'unit'
##   ..$ respect      : logi FALSE
##   ..$ rownames     : NULL
##   ..- attr(*, 'class')= chr [1:4] 'gtable' 'gTree' 'grob' 'gDesc'
##  - attr(*, 'class')= chr 'pheatmap'

重新排列行和列

mat_cluster <- mat[a$tree_row$order, a$tree_col$order]

mat_cluster

完成提取

##          sample_2   sample_3   sample_5   sample_6   sample_1    sample_4
## gene_5  2.7013864  0.5049798 -1.7796036 -0.3185864 -1.8931573 -0.13541785
## gene_2  0.3193550  0.4496518  1.7620763 -0.9312810 -1.3153020 -1.08782734
## gene_4  0.2298092 -0.3615045  1.8572989  0.7250737  1.0810986  0.70162614
## gene_1  0.3153119 -0.7730821 -0.5303812  0.5088226 -0.3286368 -0.85242874
## gene_3 -0.8220414 -1.1916559  0.2814619  1.8720241  0.6545161  0.04775437

提取聚类后的标准化矩阵

(a <- pheatmap(mat, scale='row', display_numbers = T))

直接提取不太方便。这可以自己先对数据scale标准化处理，再排序。

mat_scale <- round(t(apply(mat, 1, scale)),2)
colnames(mat_scale) <- colnames(mat)
mat_scale

最终结果

mat_cluster <- mat_scale[a$tree_row$order, a$tree_col$order]

mat_cluster

##        sample_2 sample_3 sample_5 sample_6 sample_1 sample_4
## gene_3    -0.88    -1.22     0.13     1.58     0.47    -0.08
## gene_4    -0.63    -1.42     1.53     0.03     0.50    -0.01
## gene_2     0.38     0.49     1.59    -0.67    -0.99    -0.80
## gene_1     1.04    -0.87    -0.45     1.38    -0.09    -1.01
## gene_5     1.69     0.39    -0.96    -0.10    -1.03     0.01

其他的图也都类似了，主要是获取变量的结构信息。