10xGenomics单细胞转录seurat包提取表达量画小提琴图
- 需求
- 提取表达
- 画小提琴图
- 对噪音处理
需求
其实seurat包已经有画小提琴图的函数,VlnPlot,参数也较多,基本上可以满足需求,但是如果需要对其x、y轴的文字以及大小调整,可以直接添加其他的ggplot的命名,比如命令如下:
p<-VlnPlot(object = newname.immune.combined, features.plot = as.character(gene), return.plotlist=TRUE,point.size.use = -1,cols.use=color,size.x.use = -1,do.sort = F)
p<-p+labs(x='Clusters', y= 'Gene Expresion')
p<-p+theme(panel.grid=element_blank(), legend.background = element_rect(colour = NA),
legend.title = element_blank(),legend.text = element_text(face="plain", color="black",size = 20),
axis.text.x = element_text(color="black",size=15),
axis.text.y = element_text(color="black",size=15),
axis.title.x = element_text(face="plain", color="black",size=20),
axis.title.y = element_text(face="plain", color="black",size=20))
- [x] 参数: size.x.use = -1
- [x] 参数: immune.combined为seurat对象
- [x] 参数: features.plot,为基因名称
- [x] 参数: cols.use为调整小提琴的颜色
- [x] 参数: labs添加坐标轴或者title名字,与普通的ggplot中的一样
正常情况下,这样基本上满足需求,但是之前有个特殊的项目,将项目中的聚类名字修改了,改变顺序,然后画出来的图顺序不对(如下图),因此需要调整顺序。
顺序不对的小提琴图
但是VlnPlot函数中没有找到相关参数,幸亏seurat包功能比较全,那我就提取表达量数据,自己作图就ok了。
提取表达量
seurat提取表达量函数为:FetchData,提取命令也很简单:
exprs <- data.frame(FetchData(object = immune.combined, vars.all = gene))
- [x] 参数: object,为seurat分析对象
- [x] 参数: vars.all,提取数据的名称,可以是向量,比如:c("tSNE_1","tSNE_2", 基因名称)
-
[x] 参数: cells.use,提取数据的细胞barcode,这里不用
提取表达结果如下图:
初始表达量结果
这里提取的数据只有barcode和表达量,没有聚类信息,因此需要增加聚类信息:
exprs$Barcod<-rownames(exprs)
ident<-data.frame()
#barcode与聚类信息提取
ident<-data.frame(Barcod=names(newname.immune.combined@ident),orig.ident=newname.immune.combined@ident)
#通过merge函数,将表达量与聚类号对应起来
c<-merge(exprs,ident,by='Barcod')
#对其进行排序
c$orig.ident<-factor(c$orig.ident,levels=c(sort(unique(immune.combined@ident))))
最终得到如下表达量矩阵:
有聚类信息的表达量矩阵
得到上述表达矩阵,开始画小提琴图:
ggplot(data = c,mapping = aes(x = factor(x = orig.ident),y = c[,2])) +geom_violin(scale = "width",adjust =1,trim = TRUE,mapping = aes(fill = factor(x = orig.ident)))+labs(x='Clusters', y= 'Gene Expresion',title=costmer_gene[m])+scale_colour_manual(values=mycolo)+theme(legend.position = 'none')+theme(panel.grid=element_blank(), legend.background = element_rect(colour = NA),
legend.title = element_blank(),legend.text = element_text(face="plain", color="black",size = 20),
axis.text.x = element_text(color="black",size=20),
axis.text.y = element_text(color="black",size=20),
axis.title.x = element_text(face="plain", color="black",size=20),
axis.title.y = element_text(face="plain", color="black",size=20))
ggsave('test.pdf')
错误的小提琴图
开始按照普通的画小提琴进行小提琴绘制:
ggplot(data,aes(x=factor(orig.ident),y=Cdh17,fill=orig.ident))+geom_violin(alpha=0.8,width=1)
ggsave('test.pdf')
得到如下图小提琴图,可以看到数据分别有明显从错误。
错误1小提琴图
正确的小提琴图
通过修改小提琴图画法,找到正确的小提琴图命名:
ggplot(data = data,mapping = aes(x = factor(x = orig.ident),y = Cdh17)) +geom_violin(scale = "width",adjust =1,trim = TRUE,mapping = aes(fill = factor(x = orig.ident)))
ggsave('test.pdf')
得到如下小提琴图:
小提琴图
从图中可以看出,数据分布基本正常,聚类10的表达量明显高于其他聚类,但是与我们直接用VlnPlot函数画图出来的图还是不一样,除了聚类10以外,其他表达量也较高。
最终正确小提琴图
通过查看seurat官网在GitHub上面的SingleVlnPlot源代码,然后发现其他官方对其有一个噪音处理过程,这是与VlnPlot画出来不一样的原因,然后我们也加上噪音处理,可以得到与VlnPlot一样的图。(VlnPlot函数是调用SingleVlnPlot画图)
noise <- rnorm(n = length(x = data[,c('Cdh17')])) / 100000
data[,c('Cdh17')] <- data[, c('Cdh17')] + noise
ggplot(data = data,mapping = aes(x = factor(x = orig.ident),y = Cdh17)) +geom_violin(scale = "width",adjust =1,trim = TRUE,mapping = aes(fill = factor(x = orig.ident)))
ggsave('test.pdf')
得到与VlnPlot的几乎一致的图:
最终与官方一致的小提琴图
,至于颜色,legend都可以通过ggplot2画图图层修改。
到此,画小提琴图终于完成了。
2019年3月27日
网友评论