這里我們分享一個熱圖函數(shù)树枫,是一個可視化的函數(shù)募舟,您只需要提供作圖的matrix數(shù)據(jù)即可。最近對于圖形注釋使用點比較上頭盹兢,所以我們這個函數(shù)的列注釋使用設(shè)置的·是點的注釋形狀邻梆。我們的函數(shù)解決了以下幾個問題:
數(shù)據(jù)的自動排序和列的自定義排列,數(shù)據(jù)會按照你需要的列順序自動從大到小排列绎秒,讓熱圖可視化更加美觀浦妄。
熱圖的列自動注釋,我們采用了點的注釋见芹。
legend的修改剂娄,包括刻度等修飾。
行的注釋玄呛,行的注釋需要自行 調(diào)整阅懦,函數(shù)設(shè)置了可選擇選項。
行的標(biāo)簽徘铝,在函數(shù)中我們設(shè)置了如果矩陣nrow行數(shù)大于100自動不顯示行名耳胎,應(yīng)為那樣都擠在一起也沒啥用。后續(xù)可自 行挑選需要的行進(jìn)行展示惕它。
需要注意一下怕午,函數(shù)有一個參數(shù)data_scale,假設(shè)你的數(shù)據(jù)不需要scale標(biāo)準(zhǔn)化淹魄,那么參數(shù)選擇T郁惜,作圖使用你的原始數(shù) 據(jù)。如果需要scale標(biāo)準(zhǔn)化甲锡,那么選擇F扳炬,會對你的數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化處理。
看一下函數(shù)參數(shù):(未標(biāo)注參數(shù)與Complexheatmap::Heatmap一致)
接下來我們看看函數(shù)具體的使用搔体,首先我們用一個ATAC TF分析的數(shù)據(jù),這個矩陣是已經(jīng)導(dǎo)出的半醉,行是celltype疚俱,列是TF。其 實就是我們平時做熱圖的時候一樣的數(shù)據(jù)缩多。Order_col就是設(shè)置列的順序呆奕,我這里直接演示就使用了數(shù)據(jù)的列名排列,如果 需要自己設(shè)置衬吆,這里傳入一個向量梁钾,設(shè)置順序即可。顏色也是可以設(shè)置的逊抡。
setwd("D:/KS項目/公眾號文章/scATAC-scRNA marker基因注釋熱圖")
#=================================================================================
#---------------------------------------------------------------------------------
TF <- c("JUND (264)","FOSB (190)",
"JDP2 (178)","TAL1 (42)",
"TAL2 (42)","TCF21 (21)",
"TCF23 (21)","BCL11A (799)",
"BCL11B (799)","ELF2 (364)",
"NFIX (89)","GBX2 (87)",
"SOX9 (145)","SOX4 (95)",
"SOX13 (92)","SOX12 (91)",
"NHLH2 (149)","TCF12 (159)",
"ZBTB7A (98)")
plot_atac <- read.csv('plot_atac.csv', header = T, row.names = 1)
pdf('heatmap1.pdf', width=4, height=6)
heata = ks_Heatmap(mat = plot_atac,
data_scale = T,
Order_col = colnames(plot_atac),
legendTitle = "Norm.Enrichment -log10(P-adj)[0-Max]",
point_size = 5,
heat_col = c("#E6E7E8","#3A97FF","#8816A7","black"),
TitlePosition = "leftcenter-rot",
legend_at = c(0,50,100),
legend_Lab = c(0,50,100),
column_names_gp = gpar(fontsize = 6),
row_names_gp = gpar(fontsize = 6),
customRowLabel=TF,
colanno_color = c("#7DD06F","#844081","#688EC1","#C17E73","#484125","#6CD3A7","#597873",'#3361A5'),
rowAnn = T,
rowAnno_num = c(16,3,1,4,12,9,8,11))
draw(heata, merge_legends = TRUE,heatmap_legend_side = "right")
dev.off()
接下來我們用單細(xì)胞轉(zhuǎn)錄組的數(shù)據(jù)進(jìn)行演示:我們這里的演示選擇的是每個celltype的top10marker基因姆泻。首先計算 marker基因零酪,并計算每個celltype的平均值。獲得作圖數(shù)據(jù)拇勃。
library(Seurat)
library(dplyr)
#單細(xì)胞數(shù)據(jù)的演示\先找marker基因
DefaultAssay(sce) <- "RNA"
Idents(sce) <- "celltype"
all.markers <- FindAllMarkers(sce, only.pos = TRUE, min.pct = 0.5, logfc.threshold = 0.5)
top10 <- all.markers %>% group_by(cluster) %>% top_n(n=10, wt=avg_log2FC)
table(top10$cluster)
# SMC LY UEC SF CEC EC MAC
# 10 10 10 10 10 10 10
#計算平均表達(dá)量
gene_cell_exp <- AverageExpression(sce,features = top10$gene,group.by = 'celltype',slot = 'data')
gene_cell_exp <- as.data.frame(gene_cell_exp$RNA)
作圖:
genes = c("ACTA2",
"ADIRF",
"MYL9",
"TPM2",
"CCL4",
"CCL5",
"NKG7",
"CD96",
"RHEX",
"NPAS3",
"SFRP4",
"COL3A1",
"COL1A1",
"IGF1",
"RORB",
"CAPS",
"CFAP299",
"CFAP47",
"RSPH1",
"DNAH11",
"PTPRB",
"INSR",
"ENPP2",
"IL1B",
"HLA-DQA1",
"HLA-DRA",
"HLA-DPA1",
"HLA-DRB1",
"HLA-DPB1",
"CXCL8")
#plot
pdf('heatmap2.pdf', width=4, height=6)
heata = ks_Heatmap(mat = gene_cell_exp,
data_scale = F,
Order_col = colnames(gene_cell_exp),
legendTitle = "Expression",
point_size = 5,
heat_col = c('#e0f3db','#a8ddb5','#4eb3d3','#0868ac','#084081'),
TitlePosition = "leftcenter-rot",
legend_at = c(-1,0,1,2,3),
legend_Lab = c(-1,0,1,2,3),
column_names_gp = gpar(fontsize = 6),
row_names_gp = gpar(fontsize = 6),
customRowLabel=genes,
colanno_color = c("#7DD06F","#844081","#688EC1","#C17E73","#484125","#6CD3A7","#597873"),
rowAnn = T,
rowAnno_num = c(10,10,10,10,10,10,10))
draw(heata, merge_legends = TRUE,heatmap_legend_side = "right")
dev.off()
這里我們可能會發(fā)現(xiàn)一個問題四苇,明明是每個celltype10個gene,為什么行注釋好像不是很對方咆,這是應(yīng)為有些基因不僅在一 中celltype中高表達(dá)月腋,而數(shù)據(jù)是按照表達(dá)從高到低排序的,所以才會出現(xiàn)這個問題瓣赂,可以自行調(diào)整注釋的數(shù)目榆骚,或者不采用注釋等。一個函數(shù)不可能完成所有的事情煌集,但是我們提供了框架和思路妓肢,可自行修改拓展。覺得分享有用的點個贊牙勘,分享一下再走唄职恳!