Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

这算是不一样吗 #5857

Closed
ixxmu opened this issue Oct 31, 2024 · 2 comments
Closed

这算是不一样吗 #5857

ixxmu opened this issue Oct 31, 2024 · 2 comments

Comments

@ixxmu
Copy link
Owner

ixxmu commented Oct 31, 2024

https://mp.weixin.qq.com/s/6HOyWXblLI_-EGD2KEZX_Q

@ixxmu
Copy link
Owner Author

ixxmu commented Oct 31, 2024

这算是不一样吗 by Investigator进修录

在《单细胞天地》的交流群,看到有粉丝提问关于 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数据集哦。


往期回顾

Seurat4.0系列教程11:使用sctransform

OSCA单细胞数据分析笔记13—Multi-sample comparison

scRNA-seq计算方法的优势和局限性






如果你对单细胞转录组研究感兴趣,但又不知道如何入门,也许你可以关注一下下面的课程



看完记得顺手点个“在看”哦!


生物 | 单细胞 | 转录组丨资料
每天都精彩

长按扫码可关注



1 similar comment
@ixxmu
Copy link
Owner Author

ixxmu commented Oct 31, 2024

这算是不一样吗 by Investigator进修录

在《单细胞天地》的交流群,看到有粉丝提问关于 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数据集哦。


往期回顾

Seurat4.0系列教程11:使用sctransform

OSCA单细胞数据分析笔记13—Multi-sample comparison

scRNA-seq计算方法的优势和局限性






如果你对单细胞转录组研究感兴趣,但又不知道如何入门,也许你可以关注一下下面的课程



看完记得顺手点个“在看”哦!


生物 | 单细胞 | 转录组丨资料
每天都精彩

长按扫码可关注



@ixxmu ixxmu changed the title archive_request 这算是不一样吗 Oct 31, 2024
@ixxmu ixxmu closed this as completed Oct 31, 2024
@ixxmu ixxmu closed this as completed Oct 31, 2024
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

1 participant