专栏首页生信技能树这样画热图,涉嫌操纵数据了吗

这样画热图,涉嫌操纵数据了吗

很多朋友都有这样的疑问,为什么别人绘制出来的热图,差异那么明显,除了首先他们本身就先做了差异分析,挑选出来了有差异的基因,然后才热图可视化外,其实还有一个步骤,就是按照基因(行)对表达矩阵进行zscore转换。首先看原始表达矩阵热图

代码如下:

# 2.热图
load(file='heatmap_input.Rdata')
## 2.1 数据预处理
t <- log2(cgexp+1)
t <- na.omit(t) #取出为空(表达矩阵中没有)的基因
dim(t)
## 准备画图
library(pheatmap) #加载包
identical(colnames(t),colnames(mRNA_exp_small)) #确定一下顺序没变,方便后面添加分组信息
### 构建注释矩阵
col <- data.frame(Type=group_list) #显示肿瘤类型
col$Type <- factor(col$Type,levels = c('normal','tumor'))
rownames(col) <- colnames(t)
### 开始画图

如下:

原始表达矩阵热图

可以看到,两个分组差异是有的,但是肉眼其实看不清楚基因层面哪些高表达哪些低表达。因为不同基因的表达矩阵本身差异很大,但其实我们仅仅是关心同一个基因在不同分组样本的表达,我们并不会关系不同基因的表达量问题,所以需要按照基因(行)对表达矩阵进行zscore转换。

第一次zscore

代码如下:

## 2.2 进行scale,但是不设定最大最小值
t <- t(scale(t(t)))
p1 <- pheatmap(t,
         annotation_col = col,
         annotation_legend = T,
         border_color = 'black',#设定每个格子边框的颜色,border=F则无边框
         cluster_rows = T, #对行聚类
         cluster_cols = T, #队列聚类
         show_colnames = F, #是否显示列名
         show_rownames = F #是否显示行名
)

出图如下:

第一次zscore

可以看到,上下调的基因清晰可见,但是这个时候很多人不理解热图上面的层次聚类,会错误的以为tumor被normal隔离成为了两个分组。

这个时候,如果你使用我的热图代码,通常是会在zcore的时候,设置一个上限值,比如2或者1.6,代码如下:

然后限定zscore的范围

代码如下:

## 2.3 进行scale后设定最大最小值的情况 
table(abs(t)>2)
t[t>=2]=2
t[t<=-2]=-2
table(is.na(t))
p2 <- pheatmap(t,
               annotation_col = col,
               annotation_legend = T,
               border_color = 'black',#设定每个格子边框的颜色,border=F则无边框
               cluster_rows = T, #对行聚类
               cluster_cols = T, #队列聚类
               show_colnames = F, #是否显示列名
               show_rownames = F #是否显示行名
)

出图如下:

限定zscore的范围

很有意思,这个时候上下调基因仍然是清晰可见,而且很容易看出来高低表达量分组,而且不会出现上面热图的tumor被normal隔离成为了两个分组的假象!

那么问题来了,这样是操纵数据吗

是不是很有意思,有时候你很难给合作者解释清楚。

本文分享自微信公众号 - 生信技能树(biotrainee)

原文出处及转载信息见文内详细说明,如有侵权,请联系 yunjia_community@tencent.com 删除。

原始发表时间:2020-05-17

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

我来说两句

0 条评论
登录 后参与评论

相关文章

  • 使用gtex数据库找组织特异性表达基因

    数据量可以说是很可观了,前面我们介绍过很多了:TCGA的28篇教程-GTEx数据库-TCGA数据挖掘的好帮手

    生信技能树
  • 生信分析人员如何系统入门Linux(2019更新版)

    在生信分析人员如何系统入门R(2019更新版) 里面,我提到过Linux基本上几十年都没有怎么变动过基础知识的,哪怕你现在搜索到十几年前的Linux教学视频,也...

    生信技能树
  • 一文看懂WGCNA 分析(2019更新版)

    不过,我这点战绩根本就算不上什么,其实这个WGCNA包已经是十多年前发表的了,仍然是广受好评及引用量一直在增加,破万也是指日可待。

    生信技能树
  • SpringBoot引入Bean组件的几种方式

    众所周知,在spring内部,所有的bean都是交由spring来统一管理的,有些bean可能是直接通过BeanDefinitionRegistry定义的,有些...

    奕仁
  • 小鼠大脑之空间转录组分析

    这是ST团队进行的一个示例实验。在一个单独的实验中揭示所有具有空间分辨率的小鼠大脑mRNA。

    生信技能树jimmy
  • 案例|RAC 添加表空间误将数据文件放本地处理办法

    不管是 10g、11g、12c ,RAC 架构的数据文件均是存放在共享存储上的,但是由于扩容时误操作少写一个符号,则会将数据文件建立到本地文件系统上的 $ORA...

    JiekeXu之路
  • excel分类汇总

    今天要跟分享的是excel的分类汇总功能! 分类汇总是excel中处理表格数据使用频率非常高的基础功能,可以胜任基础的统计汇总任务! 今天以一个案例作为主要介...

    数据小磨坊
  • Git和Github系列课程

    从初识git和github,到的掌握其基本使用以及使用git和github来解决多人协作问题。Git,你值得拥有。下面列举些自己收藏好的学习Git和Github...

    zhisheng
  • Git合并不同url的项目

    摘要:为了让项目能实现Git+Gerrit+Jenkin的持续集成,我们把项目从Git上迁移到了Gerrit上,发现有的同事在老Git提交代码,因为Gerrit...

    机智的程序员小熊
  • 一步一步学Python(2) 连接多台主机执行脚本

    最近在客户现场,每日都需要巡检大量主机系统的备库信息。如果一台台执行,时间浪费的就太冤枉了。 参考同事之前写的一个python脚本,配合各主机上写好的shel...

    Alfred Zhao

扫码关注云+社区

领取腾讯云代金券