这算是不一样吗

在《单细胞天地》的交流群,看到有粉丝提问关于 FindMarkers与AverageExpression 两个函数的差异 :

请问下怎么得到FindMarkers计算时ident.1和ident.2指定的cluster的平均表达量,因为输出结果中只会给avg_log2FC,我知道AverageExpression()函数可以计算平均表达量,但对应计算出来的logFC对应不上!都弄了2天了,一直找不到原因和解决办法,求大神解答!

而且还很贴心的给出来了一个示意图,如下所示:

 

这样的提问太棒了,所以我就让他把这个数据save成为rdata文件,然后把代码复制粘贴到交流群里面。我马上开始debug模式!

我的代码如下所示:

library(Seurat)
load('pbmc_Platelet_DC.Rdata')
sce=pbmc_Platelet_DC
table(sce$seurat_clusters)
table(Idents(sce))
#计算平均表达量
av <-AverageExpression(sce, assays = "RNA") 
av=as.data.frame(av$RNA)
av$diff = log(av[,2]+0.1) - log(av[,1]+0.1)

#计算差异基因
diff  <- FindMarkers(sce,
                     ident.1 = "Platelet", 
                     ident.2 = "DC",assay = "RNA")

head(diff)
head(av)

很快就能重复出来 FindMarkers与AverageExpression 两个函数的结果 :

> head(diff)
                 p_val avg_logFC pct.1 pct.2    p_val_adj
PPBP      1.511444e-10  6.028905 1.000 0.031 2.072795e-06
PF4       3.476540e-10  4.885316 1.000 0.062 4.767728e-06
GNG11     4.483201e-10  4.581688 0.929 0.000 6.148262e-06
HIST1H2AC 4.483201e-10  4.060724 0.929 0.000 6.148262e-06
SPARC     4.483201e-10  3.554182 0.929 0.000 6.148262e-06
GP9       4.483201e-10  3.553585 0.929 0.000 6.148262e-06
> head(av)
                      DC Platelet       diff
AL627309.1    0.00000000        0  0.0000000
AP006222.2    0.00000000        0  0.0000000
RP11-206L10.2 0.08462847        0 -0.6131753
RP11-206L10.9 0.00000000        0  0.0000000
LINC00115     0.00000000        0  0.0000000
NOC2L         0.32861136        0 -1.4553804

如果这样的肉眼去查看,就是粉丝提问的截图,当然是看不出来啊!会错误的以为,两个函数得到的结果差异很大!

我这里做一个小小的可视化,代码如下:


ids=intersect(rownames(diff),
              rownames(av))
comp=data.frame(FindMarkers=diff[ids,'avg_log2FC'],
           AverageExpression=av[ids,'diff'])
head(comp)
plot(comp)
head(comp)
library(ggpubr)
ggscatter(comp, x = "FindMarkers", y = "AverageExpression",
          color = "black", shape = 21, size = 3, # Points color, shape and size
          add = "reg.line",  # Add regressin line
          add.params = list(color = "blue", fill = "lightgray"), # Customize reg. line
          conf.int = TRUE, # Add confidence interval
          cor.coef = TRUE, # Add correlation coefficient. see ?stat_cor
          cor.coeff.args = list(method = "pearson",  label.sep = "\n")
)

结果如下:

 

可以看到,两个函数,FindMarkers与AverageExpression,其实计算得到 差异变化几乎是没有差异,这样的相关性已经是非常好了。(可能是粉丝对数据的一致性这个判断没有相应的统计学背景知识)

那么,有没有可能让两个函数 FindMarkers与AverageExpression得到的结果完全一致呢?

那就期待明天的讲解哈!

另外,如果你也想拿到前面的 pbmc_Platelet_DC.Rdata 文件,走这个代码,其实超级容易,就是pbmc3k数据集哦。

(0)

相关推荐

  • 单细胞工具箱|Seurat官网标准流程

    学习单细胞转录组肯定先来一遍Seurat官网的标准流程. 数据来源于Peripheral Blood Mononuclear Cells (PBMC),共2700个单细胞, Illumina Next ...

  • 工作鄙视链:公务员>事业单位>国企,私企上班只算是个打工仔

    你在北上广深的一线城市当码农,五一终于轮到你放假了,你兴奋地背上书包,穿上自己最喜欢的格子衫,坐着火车唱着歌回到了老家.你月薪2万,大城市白领,种种光环让你能够从内心俯视家乡的亲朋好友,这些土豪可能都 ...

  • 拿姿势跑法当靶子算是常态,理解

    『跑你』普及跑步/铁三/越野/急救.运动+安全=全为健康. 原创技术文章,真实案例分析 文 | 老Q 某公众号发的一篇关于跑步是否利用重力的文章,在姿势跑法教练群引发讨论. 好奇打开看了一下,其实通篇 ...

  • 再问雅与俗,诗歌算是高雅的,但废话体詩和黄诗很流行,算啥呢?

    流行就对啦,什么叫流行?流行就是人们对一种新鲜事物持有的一种短暂的狂热的追求,比如一首歌,一件衣服,一款发型等.它的特点是,短暂而汹涌,汹涌过后便被人们抛诸脑后. 与之相对的,是经典的,永恒的.不管过 ...

  • 女性若是割掉了乳房和子宫,还算是女人吗?医生给出这样的回答

    女性若是割掉了乳房和子宫,还算是女人吗?医生给出这样的回答

  • 断舍离,才算是最佳的余生

    人自呱呱落地就在不停地做着加法.学各种各样生存技能,不停累积各种各样常识,不断扩展自身的社交圈.这种过程使我们产生自身的人格,得到必要的能力,让心灵更加丰富多彩,格局更加开阔,生命更加厚实. 但伴随着 ...

  • 遇上爱情,我算是输了,经典说说

    思念是人生的一种读懂,也是一种卑微的天涯,冷风吹散,读懂是一种伤害,唯一的错,唯一的情眼,你的世界,我的清风,梦冷今生缘,彼此的伤害,只是我的读懂,我的唯一,是今生的错,也是无奈的一世梦. 今世的江湖 ...

  • 糖化血红蛋白6.5%,但空腹血糖并不高,为什么也算是糖尿病?

    有位朋友给华子留言,他去医院体检,糖化血红蛋白6.5%,但空腹血糖只有5.8mmol/L,很正常,但是医生说他是糖尿病,让他接受治疗.他有点不懂,问华子,有了糖尿病为什么空腹血糖还是正常的? 华子说, ...

  • 一碗面做成这样,算是及格了吧

    哈喽,大家好,我是吸面狂魔野食小哥. 如果这辈子只能吃一种食物,我选择面条~辣的.甜的.咸的.鲜的.麻的,兼容并包~ 今天的视频,是关于一碗不简单的面!? 碳水炸弹预警 请在食物陪伴下阅读此文 我不允 ...

  • 一碗水饺做成这样,也算是及格了

    大家好,我是野食小哥. 说起饺子,相信不论是南方人还是北方人都会对它情有独钟. 作为咱们传统美食中十分特殊的存在,饺子距今已经有1800多年的历史了. 别看它个头不大,里面却承载了幸福与团圆,让人如何 ...

  • 白酒常见的12种香型,喝过5种就算是老酒客,全喝过的不多

    我国有着几千年光辉灿烂的文明,"酒文化"也伴随着文明发展走过了数千年的历史而源远流长.我国酿酒的历史最早可追溯到新石器时代,白酒酿造以各种淀粉含量丰富的粮食为原料,以酒曲为发酵媒介 ...