最近看到一個(gè)公眾號(hào)超营,分享了一個(gè)圖片,公眾號(hào)的名字為西紅柿的空間轉(zhuǎn)錄組阅虫,代碼均來(lái)著這個(gè)公眾號(hào)演闭,廢話不多說(shuō),先看看圖長(zhǎng)什么樣
直接上代碼
加載數(shù)據(jù)與包
library(ggplot2)
library(ggrepel)
library(tidyverse)
library(Seurat)
library(SeuratData)
#載入示例數(shù)據(jù)
scRNA_harmony <- readRDS("scRNAsub.rds")
再看一下圖片颓帝,我們要構(gòu)建一個(gè)怎么樣的數(shù)據(jù) 數(shù)據(jù)里得有什么
在我看來(lái)這個(gè)圖是三個(gè)數(shù)據(jù)的體現(xiàn)
1.圖上的點(diǎn)點(diǎn) 這需要一個(gè)數(shù)據(jù) 數(shù)據(jù)1
2.圖上的中間為label的celltype需要一個(gè)數(shù)據(jù) 數(shù)據(jù)2
3.圖上那些note出來(lái)的基因需要一個(gè)數(shù)據(jù) 數(shù)據(jù)3
而他們之間為什么能相互歸屬米碰,是因?yàn)樗麄冎g有相互聯(lián)系
#開(kāi)始作圖窝革,
#先創(chuàng)建數(shù)據(jù)2 先得到celltype的label的數(shù)據(jù)
#由于我們的label是celltype 所以我們要先ident為celltype
Idents(scRNA_harmony) <- "celltype"
levels(scRNA_harmony) #查看celltype的排列順序 后面創(chuàng)建的向量要與這個(gè)排列順序?qū)?yīng)
data <- data.frame(x=0:6, #這里需要根據(jù)celltype的類型多少進(jìn)行改變
y=0,
label=c(
"Macrophages","epithelial cells","Fibroblast",
"T cells","Endothelial cells","B cells","Smooth muscle")) #數(shù)據(jù)2創(chuàng)建完成
#創(chuàng)建數(shù)據(jù)1 點(diǎn)點(diǎn)代表基因 高低代表logFC的大小 所以這個(gè)表格大致就是FindAllMarkers生成的表格
marker<-FindAllMarkers(scRNA_harmony)
#再觀察圖片,發(fā)現(xiàn)點(diǎn)點(diǎn)(基因)是隨機(jī)排列的见间,所以我們要給每一個(gè)基因弄一個(gè)隨機(jī)數(shù)聊闯,好讓他們隨機(jī)排列
#第一步 在marker這個(gè)表格中給每一種celltype按照排序給他們定義一個(gè)數(shù)字0-....
for(i in 1:nrow(data)){
marker[which(marker$cluster == data$label[i]),'id'] <- data$x[i]
}
marker$id=as.numeric(marker$id)
#第二步 #根據(jù)這個(gè)數(shù)字來(lái)創(chuàng)建隨機(jī)數(shù)
marker$x=0
for (i in 0:9) {
marker[which(marker$id == i),"x"]= i+runif(nrow(marker[which(marker$id == i),]),
min = -0.4, #為啥0.4,因?yàn)?/2=0.5
max = 0.4)
} #這樣數(shù)據(jù)1就創(chuàng)建好了
#創(chuàng)建數(shù)據(jù)3 創(chuàng)建這個(gè)表格的目的就是顯示top基因的名字
#所以是基因marker這個(gè)數(shù)據(jù) 來(lái)得到top基因的
marker$p=ifelse(marker$p_val_adj<0.01,"a","b")#設(shè)置一個(gè)p值的顏色label
top = marker %>% group_by(cluster) %>% top_n(n = 5, wt = avg_log2FC)#前5
bt= marker %>% group_by(cluster) %>% top_n(n = 5, wt = -avg_log2FC)#后5
xx=rbind(top,bt) #這樣數(shù)據(jù)3就創(chuàng)建好了
#數(shù)據(jù)1-marker
#數(shù)據(jù)2-data
#數(shù)據(jù)3-xx
數(shù)據(jù)已經(jīng)創(chuàng)建好了 現(xiàn)在開(kāi)始畫圖
#首先米诉,設(shè)置顏色菱蔬,該顏色的排列順序就是下面label celltype的顏色 所以可以把這個(gè)排列順序調(diào)成跟seurat中的celltype的顏色一樣
mycol<-c("#223D6C","#D20A13","#FFD121","#088247",
"#11AA4D","#58CDD9","#7A142C","#5D90BA",
"#431A3D","#91612D","#6E568C","#E0367A",
"#D8D155","#64495D","#7CC767")
#畫圖 先畫數(shù)據(jù)1-marker 主圖
p<-ggplot(data = marker, mapping = aes(x = x,
y = avg_log2FC)) +
geom_point(aes(colour = p),size=0.4) +
scale_color_manual(values = c("red","black"),
labels = c(paste('Adjusted P-val <',0.01),
paste('Adjusted P-val >=',0.01)))+
xlab("") + ylab("average logFC") +
geom_tile(data = data.frame(x=(0:6),
y=0,
label=(0:6)),
aes(x=x,y=y,fill=factor(label),
height = 0.4))+
scale_fill_manual(values = mycol)
p
#加上數(shù)據(jù)2的label
p2=p+geom_text(data=data ,
aes(x=x,y=y,label=label),
color="white")+ #改細(xì)胞名字體的顏色 想一下該字體的大小呢 #足以見(jiàn)系統(tǒng)學(xué)習(xí)gglot有多么重要
geom_text_repel(data=xx,
aes(x = x,
y = avg_log2FC,
label=gene),
size=3)
p2
#加上數(shù)據(jù)3的top基因
p2+theme_classic()+
guides(fill = FALSE,size=FALSE)+
theme(legend.position = "top",
axis.text.x = element_blank(),
axis.line.x = element_blank(),
axis.ticks.x = element_blank())
馬上要放五一了,如坐針氈史侣,心早就飛了拴泌,所以來(lái)寫兩篇文章消耗一下時(shí)間,放假前的打算(系統(tǒng)的學(xué)一下ggplot)惊橱,我覺(jué)得是不可能的........我........
References:跟著Cell學(xué)畫圖~ (qq.com)