我们之前写过计算单细胞比例的很多帖子,目的就是看看对照组与实验组相比,是否有些细胞类型或者亚群数量发生了显著改变:
最近,看到之前张泽明教授团队发表在Science上的文章:《Pan-cancer single cell landscape of tumor-infiltrating T cells》,提出了一种方法,计算OR值,用来评估单细胞的组织偏好性。张教授是单细胞领域的大佬,发表在science上的方法可放心使用。此外,这篇文章的所有分析代码都是公开的,但是如果要学习还是有一定难度的,尤其是根据自己数据。
OR的方法在原文方法中描述了,这里我摘录下来,并用我的塑料英语水平进行了翻译,应该没有问题,帮助大家理解。
原文Method:To characterize the tissue distribution of meta-clusters, odds ratios (OR) were calculated and used to indicate preferences. Specifically, for each combination of meta-cluster i and tissue j, a 2 by 2 contingency table was constructed, which contained the number of cells of meta-cluster i in tissue j, the number of cells of meta-cluster i in other tissues, the number of cells of non-i meta-clusters in tissue j, the number of cells of non-i meta-clusters in other tissues. Then Fisher’s exact test was applied on this contingency table, thus OR and corresponding p-value could be obtained. P-values were adjusted using the BH method implemented in the R function p.adjust. We found that all ORs > 1.5 or ORs < 0.5 had adjusted p-values < 1e-10. Hence, a higher OR with a value > 1.5 indicated that meta-cluster i was more preferred to distribute in tissue j, a lower OR with a value < 0.5 indicated that meta-cluster i was preferred not to distribute in tissue j.
译文方法:为了表征meta-clusters的组织分布(这里的meta-cluster就是celltype或者sub-celltype),计算odds ratios(OR)比值用于指示偏好性(或者说倾向分布)。具体而言,对于每一个meta-cluster(celltype) i 和其对应的组织j,构建一个2 × 2列联表。此表包含meta-cluster(celltype) i 在组织j中、以及其他组织中的细胞数量。组织j中非meta-cluster(celltype) i的细胞数,其他组织中非meta-cluster(celltype) i的细胞数。然后对该列联表进行Fisher精确检验,得到OR和相应的p值。采用R中p.adjust函数对p值进行BH矫正。我们发现所有ORs > 1.5或ORs < 0.5的校正p值均< 1e-10。因此,较高的OR值> 1.5表示meta-cluster(celltype) i 更倾向于分布在组织j中,较低的OR值< 0.5表示meta-cluster(celltype) i 不倾向于分布在组织j中。
这个分析涉及两个函数,一个是分析函数,另一个是分析及作图函数。我们做的事情是将这个函数进行了注释,并对函数进行了改造,适用于所有情况,原文只适用于自己的分组和相应的列。我们改造后里面的参数可随意变化使用,可以使用不同的分组。
接下来看看具体的分析,加载函数和数据,对了,以下这些包都要有!
library("sscVis")
library("data.table")
library("grid")
library("cowplot")
library("ggrepel")
library("readr")
library("plyr")
library("ggpubr")
library("ggplot2")
library("dplyr")
library("tidyr")
setwd("D:/KS项目/公众号文章/张泽明code学习")
#加载函数及数据
mouse_data <- readRDS("./mouse_data.rds")
source("./test_function.R")
source("./draw_analysis.R")
数据分析并做图,分析文件存储在list中。
#数据分析
A <- do.tissueDist(cellInfo.tb = mouse_data@meta.data,#这里的输入文件需要的是整理好的有分组和细胞类型的metadata文件
out.prefix = "./ORplot",#设定导出的文件名,自己设置
pdf.width = 5,#热图的宽设置
pdf.height = 8,#热图的高度设置
verbose=1,#设置为1文件以list形式存储
meta.cluster = 'celltype',#这里是细胞类型,也可以是seurat_clusters,名称没有要求,就是你细胞类型的列名
loc = 'SingleCellNet_Xie',#这里就是分组,metadata中分组的列名,至于命名没有要求
z.hi=10)#热图legend最大值,根据实际情况自己设置
#查看并保存文件
A$OR.dist.mtx #做热图数据,OR值
A$p.dist.tb #p值
A$OR.dist.tb #OR值
A$count.dist.melt.ext.tb#组合表,adjust-pvalue等
分析完成拿到数据,至于作图嘛,完全可以不用函数中提供的默认作图(不过默认的图也很好),拿到矩阵后用ggplot2作图完全没有问题!由于我的示例数据很随意,不太好,所以结果一般。想要原文的效果,可以试试自己的数据或者使用原文数据!
#自己做图
data <- A$count.dist.melt.ext.tb
write.csv(data, file = 'data.csv')
library(ggplot2)
library(RColorBrewer)
data <- read.csv("data.csv", header = T)
ggplot(data, aes(cid, rid)) +
geom_tile(aes(fill = OR), colour = "black", size = 0.6)+
scale_fill_gradientn(name='OR',
colours=colorRampPalette(rev(brewer.pal(n = 7, name ="RdYlBu")))(100))+
theme_minimal() +
theme(axis.title.x=element_blank(),
axis.ticks.x=element_blank(),
axis.title.y=element_blank(),
axis.text.y = element_text(size = 10,color = 'black')) +
scale_y_discrete(position = "right")
以上就是今天所有内容了,觉得有用的点个赞再走呗!!!
更多精彩内容请至我的公众号---KS科研分享与服务