10xGenomics單細胞轉(zhuǎn)錄seurat包提取表達量畫小提琴圖
- 需求
- 提取表達
- 畫小提琴圖
- 對噪音處理
需求
其實seurat包已經(jīng)有畫小提琴圖的函數(shù)恬惯，VlnPlot，參數(shù)也較多，基本上可以滿足需求家破，但是如果需要對其x韭寸、y軸的文字以及大小調(diào)整，可以直接添加其他的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] 參數(shù)： size.x.use = -1
- [x] 參數(shù)： immune.combined為seurat對象
- [x] 參數(shù)： features.plot趟大，為基因名稱
- [x] 參數(shù)： cols.use為調(diào)整小提琴的顏色
- [x] 參數(shù)： labs添加坐標軸或者title名字，與普通的ggplot中的一樣
正常情況下铣焊，這樣基本上滿足需求逊朽，但是之前有個特殊的項目，將項目中的聚類名字修改了曲伊，改變順序惋耙，然后畫出來的圖順序不對（如下圖），因此需要調(diào)整順序熊昌。
順序不對的小提琴圖
但是VlnPlot函數(shù)中沒有找到相關(guān)參數(shù)绽榛，幸虧seurat包功能比較全，那我就提取表達量數(shù)據(jù)婿屹，自己作圖就ok了灭美。
提取表達量
seurat提取表達量函數(shù)為：FetchData，提取命令也很簡單：
exprs <- data.frame(FetchData(object = immune.combined, vars.all = gene))
- [x] 參數(shù)： object昂利，為seurat分析對象
- [x] 參數(shù)： vars.all届腐，提取數(shù)據(jù)的名稱铁坎，可以是向量，比如：c("tSNE_1","tSNE_2", 基因名稱)
-
[x] 參數(shù)： cells.use犁苏，提取數(shù)據(jù)的細胞barcode硬萍，這里不用
提取表達結(jié)果如下圖：
初始表達量結(jié)果
這里提取的數(shù)據(jù)只有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函數(shù)朴乖，將表達量與聚類號對應(yīng)起來
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')
得到如下圖小提琴圖助赞，可以看到數(shù)據(jù)分別有明顯從錯誤买羞。
正確的小提琴圖
通過修改小提琴圖畫法，找到正確的小提琴圖命名：
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')
得到如下小提琴圖：
從圖中可以看出雹食，數(shù)據(jù)分布基本正常畜普，聚類10的表達量明顯高于其他聚類，但是與我們直接用VlnPlot函數(shù)畫圖出來的圖還是不一樣群叶，除了聚類10以外吃挑，其他表達量也較高。
最終正確小提琴圖
通過查看seurat官網(wǎng)在GitHub上面的SingleVlnPlot源代碼街立，然后發(fā)現(xiàn)其他官方對其有一個噪音處理過程儒鹿，這是與VlnPlot畫出來不一樣的原因，然后我們也加上噪音處理几晤，可以得到與VlnPlot一樣的圖约炎。（VlnPlot函數(shù)是調(diào)用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日
