隨著時(shí)代的發(fā)展,單一研究轉(zhuǎn)錄組、蛋白代謝、甲基化等已經(jīng)難以滿足研究者越來越高的研究期望,大家更多地期望聯(lián)合多種數(shù)據(jù)進(jìn)行多組學(xué)聯(lián)合分析。那么這時(shí)候,一種好的展示結(jié)果的方式無疑會(huì)為發(fā)表高分文章增光添彩。
誠然,小編的上一篇(如何讓你的圖變得高大上之ComplexHeatmap())結(jié)局得太匆匆,因此,本篇奉上作圖代碼實(shí)例。
本次,我們將展示一個(gè)甲基化與表達(dá)譜聯(lián)合分析的熱圖。本著先學(xué)習(xí)再創(chuàng)造的態(tài)度,小編做了一下知識(shí)的搬運(yùn)工,本篇所有代碼均引用自:Zuguang Gu, Roland Eils and Matthias Schlesner,?Complex heatmaps reveal patterns and correlations in multidimensional genomic data,?Bioinformatics, 2016。
我們先上效果圖:
其實(shí)代碼并不長(zhǎng),關(guān)鍵在于如何準(zhǔn)備作圖所用數(shù)據(jù)以及對(duì)代碼的理解上,所以下面將著重對(duì)這兩點(diǎn)進(jìn)行解釋說明。
>library(ComplexHeatmap)
>library(circlize)
#為了測(cè)試方便,用setwd設(shè)置工作路徑
setwd(“/media/bmk/***/Meth_Exp”)
#rds格式是一種保存了數(shù)據(jù)的數(shù)據(jù)類型和結(jié)構(gòu)的格式,可以用函數(shù)saveRDS()來生成這種格式的文件。meth.rds文件來自R包示例文件,如需查看這個(gè)文件中都是什么,請(qǐng)參考(看不清熱圖的聚類結(jié)果怎么辦)的方法。另外,文件中只包含了DMR關(guān)聯(lián)基因的甲基化和表達(dá)量呈現(xiàn)負(fù)相關(guān)的DMR
>res_list?=?readRDS(“data/meth.rds”)
#顯示樣本是腫瘤樣本還是正常樣本的標(biāo)簽
>type?=?res_list$type
(20個(gè)正常樣本和20個(gè)腫瘤樣本)
#矩陣,其中行對(duì)應(yīng)差異甲基化區(qū)域(DMR),矩陣中的值是每個(gè)樣品中DMR中的平均甲基化水平。
>mat_meth?=?res_list$mat_meth
#矩陣,其中行對(duì)應(yīng)于與DMR相關(guān)的基因(即與DMR最近的基因,也可以叫做DMR關(guān)聯(lián)基因),矩陣中的值是每個(gè)樣品中每個(gè)基因的表達(dá)水平(對(duì)樣品中的每個(gè)基因的表達(dá)進(jìn)行了標(biāo)準(zhǔn)化)。
>mat_expr?=?res_list$mat_expr
#甲基化變化的方向(hyper表示腫瘤樣本中高甲基化,hypo表示腫瘤樣本中低甲基化)
>direction?=?res_list$direction
#甲基化與甲基化關(guān)聯(lián)基因相關(guān)性p值
>cor_pvalue?=?res_list$cor_pvalue
#基因類型(如蛋白編碼基因或lincRNA等)
>gene_type?=?res_list$gene_type
#DMR注釋到基因的功能區(qū)間(如intergenic/intragenic或者TSS
>anno_gene?=?res_list$anno_gene
#DMR到關(guān)聯(lián)基因TSS的距離
>dist?=?res_list$dist
#與增強(qiáng)子重疊的DMR的部分
>anno_enhancer?=?res_list$anno_enhancer
#顏色定義參見上一篇(如何讓你的圖變得高大上之ComplexHeatmap()) ##我們首先定義兩個(gè)列注釋,然后制作復(fù)雜的熱圖。 >ht_global_opt( #Heatmap()實(shí)際上是單一熱圖的類構(gòu)造函數(shù)。如果需要組合超過一個(gè)熱圖,用戶可以通過+操作符添加熱圖。默認(rèn)情況下,將兩個(gè)熱圖通過+連接后,第二個(gè)熱圖的行聚類樹會(huì)去掉,行的順序會(huì)與是第一個(gè)熱圖的順序保持一致。
##首先計(jì)算甲基化矩陣的列聚類,以便可以將表達(dá)矩陣中的列調(diào)整為具有與甲基化矩陣中相同的列順序。
>column_tree?=?hclust(dist(t(mat_meth)))
>column_order?=?column_tree$order
>library(RColorBrewer)
#定義甲基化表達(dá)水平顏色,從0/blue-0.5/white-1/red漸變
>meth_col_fun?=?colorRamp2(c(0,?0.5,?1),?c(“blue”,?“white”,?“red”))
#定義甲基化變化方向?qū)?yīng)顏色
>direction_col?=?c(“hyper”?=?“red”,?“hypo”?=?“blue”)
#定義表達(dá)水平顏色
>expr_col_fun?=?colorRamp2(c(-2,?0,?2),?c(“green”,?“white”,?“red”))
#定義相關(guān)性p值顏色
>pvalue_col_fun?=?colorRamp2(c(0,?2,?4),?c(“white”,?“white”,?“red”))
#定義基因類型顏色
>gene_type_col?=?structure(brewer.pal(length(unique(gene_type)),?“Set3”),
names?=?unique(gene_type))
#定義注釋model顏色
>anno_gene_col?=?structure(brewer.pal(length(unique(anno_gene)),?“Set1”),
names?=?unique(anno_gene))
#定義距離顏色
>dist_col_fun?=?colorRamp2(c(0,?10000),?c(“black”,?“white”))
#定義增強(qiáng)子相關(guān)顏色
>enhancer_col_fun?=?colorRamp2(c(0,?1),?c(“white”,?“orange”))
#ht_global_opt()是一個(gè)可選函數(shù),它會(huì)全局控制一些參數(shù)。我們可以通過此全局函數(shù)同時(shí)為所有熱圖/注釋設(shè)置一些參數(shù)。需要注意的是,一定將它放在熱圖代碼(也就是Heatmap())之前,并在繪制熱圖后重置所有選項(xiàng)值以消除對(duì)下一個(gè)熱圖的影響。
#可以通過?ComplexHeatmap::ht_global_opt查看此函數(shù)的幫助
>?names(ht_global_opt())????????#可查看該函數(shù)可定義的參數(shù)
[1]?“heatmap_row_names_gp”
[2]?“heatmap_column_names_gp”
[3]?“heatmap_row_title_gp”
[4]?“heatmap_column_title_gp”
[5]?“heatmap_legend_title_gp”
[6]?“heatmap_legend_title_position”
[7]?“heatmap_legend_labels_gp”
[8]?“heatmap_legend_grid_height”
[9]?“heatmap_legend_grid_width”
[10]?“heatmap_legend_grid_border”
[11]?“annotation_legend_title_gp”
[12]?“annotation_legend_title_position”
[13]?“annotation_legend_labels_gp”
[14]?“annotation_legend_grid_height”
[15]?“annotation_legend_grid_width”????[16]?“annotation_legend_grid_border”
[17]?“fast_hclust”
heatmap_legend_title_gp?=?gpar(fontsize?=?8,?fontface?=?“bold”),
heatmap_legend_labels_gp?=?gpar(fontsize?=?8),
heatmap_column_names_gp?=?gpar(fontsize?=?8),
heatmap_column_title_gp?=?gpar(fontsize?=?10),
heatmap_row_title_gp?=?gpar(fontsize?=?8)
)
#利用HeatmapAnnotation()對(duì)行或列注釋。HeatmapAnnotation()函數(shù)會(huì)生成一個(gè)注釋用的list對(duì)象。該函數(shù)的主要格式是HeatmapAnnotation(df/數(shù)據(jù)框, name/注釋名稱, col/注釋顏色列表, show_legend/是否顯示數(shù)據(jù)框中每一列的圖例)
#樣本類型注釋,Tumor樣本為pink,Control樣本為royalbule,名稱在左邊
>ha?=?HeatmapAnnotation(type?=?type,
col?=?list(type?=?c(“Tumor”?=?“pink”,?“Control”?=?“royalblue”)),
annotation_name_side?=?“left”)
#不顯示圖例名稱
>ha2?=?HeatmapAnnotation(type?=?type,
col?=?list(type?=?c(“Tumor”?=?“pink”,?“Control”?=?“royalblue”)),
show_legend?=?FALSE)
>ht_list?=?Heatmap(mat_meth,?name?=?“methylation”,?col?=?meth_col_fun,
column_order=?column_order,
top_annotation?=?ha,?column_title?=?“Methylation”)?+
#方向部分
Heatmap(direction,?name?=?“direction”,?col?=?direction_col)?+
#表達(dá)部分
Heatmap(mat_expr[,?column_tree$order],?name?=?“expression”,
col?=?expr_col_fun,
column_order?=?column_order,
top_annotation?=?ha2,?column_title?=?“Expression”)?+
#p值部分
Heatmap(cor_pvalue,?name?=?“-log10(cor_p)”,?col?=?pvalue_col_fun)+
#基因類型部分
Heatmap(gene_type,?name?=?“gene?type”,?col?=?gene_type_col)?+
#基因注釋部分
Heatmap(anno_gene,?name?=?“anno_gene”,?col?=?anno_gene_col)?+
#距離部分
Heatmap(dist,?name?=?“dist_tss”,?col?=?dist_col_fun)?+
#增強(qiáng)子部分
Heatmap(anno_enhancer,?name?=?“anno_enhancer”,?col?=?enhancer_col_fun,
cluster_columns?=?FALSE,?column_title?=?“Enhancer”)
#以上熱圖相加后得到的是
#可以看到,整個(gè)熱圖是單個(gè)熱圖從左到右依次相加,圖例是從上到下依次相加,除了第一個(gè)熱圖的行聚類樹保留了之外,其他的都默認(rèn)被去除掉,順序和第一個(gè)熱圖保持一致。
#熱圖相加的返回值是一個(gè)HeatmapList對(duì)象。直接允許ht_list對(duì)象會(huì)默認(rèn)調(diào)用draw()方法。通過顯式地調(diào)用draw()方法,你可以進(jìn)行更多的控制,例如圖例和標(biāo)題。
#可以通過??`draw,HeatmapList-method`來查看draw的幫助
#分兩個(gè)cluster進(jìn)行聚類,按hyper和hypo進(jìn)行分裂,加上整張熱圖的名稱,設(shè)置圖例位置等
>draw(ht_list,?km?=?2,?split?=?direction,
column_title?=?“Comprehensive?correspondence?between?methylation,
expression?and?other?genomic?features”,
column_title_gp?=?gpar(fontsize?=?12,?fontface?=?“bold”),
merge_legends?=?TRUE,?heatmap_legend_side?=?“bottom”)
#重置全局參數(shù)消除影響
>ht_global_opt(RESET?=?TRUE)
復(fù)雜的熱圖顯示高度甲基化的DMR富含基因間和基因內(nèi)區(qū)域,很少與增強(qiáng)子重疊。相反,低甲基化的DMR富含轉(zhuǎn)錄起始位點(diǎn)(TSS)和增強(qiáng)子。
知識(shí)點(diǎn)總結(jié)
1.ComplexHeatmap可實(shí)現(xiàn)單個(gè)熱圖的相加以實(shí)現(xiàn)數(shù)據(jù)之間的聯(lián)合。
2.ht_global_opt()函數(shù)可實(shí)現(xiàn)整個(gè)熱圖的全局控制,但要注意使用結(jié)束后進(jìn)行重置。
3.draw()函數(shù)在返回值是HeatmapList對(duì)象可以實(shí)現(xiàn)更多的控制。