一、熱圖的介紹
(1)什么是熱圖
熱圖是heatmap的直譯揽祥,用暖色表示數(shù)值大讽膏,冷色表示數(shù)值小。行為基因拄丰,列為樣本府树。
(2)熱圖的作用
- 數(shù)據(jù)質(zhì)量控制:通過熱圖,我們能很直觀的看到不同分組之間的基因整體表達模式料按,因此可以迅速地判斷同組之間各樣本的重復(fù)性如何奄侠,從而判斷實驗處理是否正確,數(shù)據(jù)是否可信可靠符合邏輯载矿。相同實驗分組的幾個生物學(xué)重復(fù)應(yīng)該聚類在一起垄潮。
- 直觀展示重點研究對象的表達量數(shù)據(jù)差異變化情況:一次實驗中檢測到的基因往往成千上萬,導(dǎo)致一幅全局性熱圖的行數(shù)也相應(yīng)地十分龐大,單個基因的表達情況幾乎無法辨識弯洗,所以用整個數(shù)據(jù)集畫出的熱圖往往只能用于數(shù)據(jù)的整體質(zhì)控旅急。而要能夠向讀者清晰展示自己所研究的某一批基因的數(shù)據(jù)分布與變化情況,就應(yīng)當把自己的研究對象從數(shù)據(jù)集中找出來再繪制熱圖牡整,作為重點關(guān)注的幾個標志性的基因在右側(cè)標注藐吮。
- 不要畫“沒用”的熱圖,比如圖形的展示結(jié)果是早可預(yù)見的逃贝,沒用什么意義的熱圖炎码;用于進行探索性數(shù)據(jù)分析、提供下游分析處理線索的熱圖秋泳。應(yīng)該能夠解釋:這個熱圖展示的特定的表達模式潦闲,說明了哪些生物學(xué)功能與信號通路的改變?是什么導(dǎo)致了這些改變迫皱?這種改變又會產(chǎn)生什么樣的后果歉闰?
二、熱圖中數(shù)據(jù)的標準化
(1)為什么要標準化
熱圖的圖例是以0為中心卓起,變化范圍在±3左右和敬。RNA-seq得到的基因表達量其實差異是非常大的,從個位數(shù)到上千戏阅,所以首先要對基因的表達量進行歸一化昼弟,不然我們得需要多大尺度的色階才能表示這么大范圍的變化啊。標準化可以使數(shù)據(jù)范圍縮小奕筐,便于用合適的顏色色度進行表示舱痘,能避免由于單個數(shù)據(jù)過大(過小)离赫,導(dǎo)致冷暖色分布不明顯的現(xiàn)象芭逝。
(2)怎樣標準化
可以按基因、按樣本或全局標準化渊胸。一般選擇按行(基因)進行標準化旬盯。如果對同一行進行標準化,就不能比較同一樣品中不同基因的表達(不能比較不在同一行的z-score數(shù)值翎猛,原來不同基因絕對表達量高的也有可能在同一樣品中得到低的z-score分數(shù))胖翰。如果沒有異常值,也可以進行全局的標準化切厘。
(3)標準化方法
-
Z-score scaling(正態(tài)標準化):它代表原始數(shù)值和總體均值之間的距離萨咳,并以標準差為單位計算。在原始數(shù)值低于平均值時Z則為負數(shù)迂卢,反之則為正數(shù)某弦,是一個不受原始測量單位影響的數(shù)值桐汤,經(jīng)常用在數(shù)據(jù)分化比較大的場景。在分類靶壮、聚類怔毛、PCA算法中,使用z-score值的結(jié)果更好腾降。轉(zhuǎn)換后的數(shù)據(jù)符合正態(tài)分布拣度,平均數(shù)為0,方差為1螃壤。Z-score要求原始數(shù)據(jù)的分布為正態(tài)分布或者近似正態(tài)分布抗果。
z-score = (x-μ)/SD
- log對數(shù)變換:原始的counts矩陣是一個離散型的變量,離散程度很高奸晴。即使經(jīng)過了RPKM/FPKM等方法抵消了一些測序技術(shù)誤差的影響冤馏,但高低豐度基因的表達量的差距依然很大。常用的歸一化方法是對數(shù)轉(zhuǎn)換log2(counts+1)寄啼。加1是因為表達量為0的無法取到對數(shù)逮光,而加1后log2(1)還是等于0。轉(zhuǎn)換后不改變數(shù)據(jù)原本的排序關(guān)系墩划,可以在全局進行比較涕刚。
- min-max標準化(0-1標準化):是對原始數(shù)據(jù)的線性變換,使結(jié)果落到[0 , 1]區(qū)間乙帮。更符合直觀認知噪猾,使同一行里基因表達最高值為1捧搞,最低值為0处坪。
三锅必、熱圖中數(shù)據(jù)的聚類
(1)什么是聚類
聚類,就是把相似的東西分到一組塞绿。一個聚類算法通常只需要知道如何計算相似度就可以開始工作了沟涨,是無監(jiān)督學(xué)習(xí)。聚類本質(zhì)上是集合劃分問題异吻。因為沒有人工定義的類別標準,因此算法要解決的核心問題是如何定義簇喜庞,唯一的要求是簇內(nèi)的樣本盡可能相似诀浪。通常的做法是根據(jù)簇內(nèi)樣本之間的距離,或是樣本點在數(shù)據(jù)空間中的密度來確定延都。對簇的不同定義可以得到各種不同的聚類算法雷猪。常見的聚類算法有:
- 連通性聚類。典型代表是層次聚類算法晰房,它根據(jù)樣本之間的聯(lián)通性來構(gòu)造簇求摇,所有聯(lián)通的樣本屬于同一個簇射沟。
- 基于質(zhì)心的聚類。典型代表是k均值算法与境,它用一個中心向量來表示這個簇验夯,樣本所屬的簇由它到每個簇的中心距離確定。
(2)層次聚類
找出哪兩個基因最相似(similar)摔刁,將它們合并作為一個整體(merge them to a cluster)挥转;再找出哪兩個基因最相似,合并成一個整體共屈;重復(fù)這個步驟绑谣,直到所有的基因都聚成一類。
層次聚類通常會伴隨著產(chǎn)生一個系統(tǒng)樹圖(dendrogram)拗引,表示了相似性和聚類的順序借宵。
(3)相似性如何定義
可以有多種計算距離的方式,如Euclidean distance(歐氏距離)矾削、Manhattan distance等暇务,現(xiàn)實中絕大多數(shù)時候會使用歐氏距離。使用包含歐氏距離的聚類算法時怔软,要注意數(shù)據(jù)的量綱影響垦细,有必要時先進行數(shù)據(jù)標準化處理。
You need to make sure the scales for each variable are roughly equivalent, otherwise you will be biased towards one of them. The standard practice is to divide each variable by its standard deviation.
(4)聚類如何作為一個整體
- the average(類平均法):類與另一個類中任兩個樣品距離的平均值挡逼。效果比較好括改。
- the closest point:single-linkage
- the furthest point:complete-linkage
- centroid(重心):類的均值與另一個類的均值的距離
- median(中間距離法)
- ward(離差平方和法)
- density(密度估計法)
四、pheatmap(Pretty Heatmaps)
(1)R Script
> library(pheatmap)
> library(RColorBrewer)
> #測試數(shù)據(jù)
> testdata1[1:5,1:5]
CF1 CF2 CF3 CF4 CM1
rna13468 66.97984 80.07318 54.87525 91.65463 1.401584
rna885 26.53467 44.51450 33.65076 40.72113 60.633389
rna32332 0.00000 37.71825 11.89813 0.00000 1.403419
rna8744 42.44415 52.31791 60.35968 54.00533 39.524090
rna16488 0.00000 0.00000 0.00000 0.00000 0.000000
> #樣本信息
> anno_col <- data.frame(Group=c(rep("CF",4),rep("CM",4),
+ rep("PF",4),rep("PM",4)),
+ row.names = colnames(testdata1))
> anno_col$Group <- factor(anno_col$Group,levels=unique(anno_col$Group))
> #基因注釋
> head(anno_row)
gene_name chromosome
rna10784 CG13244 2L
rna12885 CG14589 2R
rna13468 CG1882 2R
rna16488 CG8311 2R
rna18953 CG5539 2R
rna20313 CG1317 3L
> anno_row$chromosome <- droplevels(anno_row$chromosome)
> #注釋顏色
> Group_col <- brewer.pal(8,"Set2")[2:5]
> names(Group_col) <- levels(anno_col$Group)
> Chr_col <- brewer.pal(8,"Set3")[1:5]
> names(Chr_col) <- levels(anno_row$chromosome)
> ann_colors <- list(Group=Group_col,
+ Chromosome=Chr_col)
> #默認參數(shù)繪圖
> pheatmap(testdata1)
> #添加注釋等
pheatmap( testdata1,
+ scale="row",
+ clustering_method="average",
+ color=colorRampPalette(c("green", "black","red"))(100),
+ annotation_col=anno_col,
+ annotation_row=anno_row[2],
+ annotation_names_row=F, annotation_names_col=T,
+ annotation_colors=ann_colors,
+ border_color = NA)
(2)用法
pheatmap(mat, color = colorRampPalette(rev(brewer.pal(n = 7, name = "RdYlBu")))(100), kmeans_k = NA, breaks = NA, border_color = "grey60", cellwidth = NA, cellheight = NA, scale = "none", cluster_rows = TRUE, cluster_cols = TRUE, clustering_distance_rows = "euclidean", clustering_distance_cols = "euclidean", clustering_method = "complete", clustering_callback = identity2, cutree_rows = NA, cutree_cols = NA, treeheight_row = ifelse((class(cluster_rows) == "hclust") || cluster_rows, 50, 0), treeheight_col = ifelse((class(cluster_cols) == "hclust") || cluster_cols, 50, 0), legend = TRUE, legend_breaks = NA, legend_labels = NA, annotation_row = NA, annotation_col = NA, annotation = NA, annotation_colors = NA, annotation_legend = TRUE, annotation_names_row = TRUE, annotation_names_col = TRUE, drop_levels = TRUE, show_rownames = T, show_colnames = T, main = NA, fontsize = 10, fontsize_row = fontsize, fontsize_col = fontsize, angle_col = c("270", "0", "45", "90", "315"), display_numbers = F, number_format = "%.2f", number_color = "grey30", fontsize_number = 0.8 * fontsize, gaps_row = NULL, gaps_col = NULL, labels_row = NULL, labels_col = NULL, filename = NA, width = NA, height = NA, silent = FALSE, na_col = "#DDDDDD", ...)
(3)參數(shù)
- scale="row" #標準化z-score,按行(基因)row/按列column/不標準化none
- cluster_rows=T, cluster_cols=T #默認聚類
- clustering_distance_rows="euclidean", clustering_distance_cols="euclidean" #默認歐氏距離
- clustering_method="complete" #默認最遠距離
- color=colorRampPalette(c("green", "black","red"))(100) #紅黑綠
- legend=TRUE, legend_breaks=NA, legend_labels=NA #默認圖例
- annotation_col=annotation_col, annotation_row=annotation_row, #行和列的注釋
- annotation_colors=ann_colors #設(shè)置分組顏色
- annotation_legend=T #注釋的圖例
- annotation_names_row=T, annotation_names_col=F #注釋的名稱
- show_rownames=T, show_colnames=T #是否顯示行名和列名
- border_color = NA #方格邊框顏色,默認"grey60" ,NA為不畫邊框
- treeheight_row=30, treeheight_col=20 #進化樹高度,設(shè)為0可隱藏進化樹
五家坎、ComplexHeatmap
(1)R Script
> library(ComplexHeatmap)
> library(circlize)
> library(RColorBrewer)
> #提前標準化數(shù)據(jù)
> heatmapdata <- testdata1
> heatmapdata <- t(scale(t(heatmapdata)))
> #heatmapdata[heatmapdata>2]=2
> #heatmapdata[heatmapdata<-2]=-2
默認參數(shù)繪圖
> Heatmap(testdata1)
> Heatmap(heatmapdata)
熱圖顏色
> col_fun = colorRamp2(c(-2, 0, 2), c("green", "white", "red"))
> Heatmap(heatmapdata, name=" ", col=col_fun)
> Heatmap(heatmapdata, name=" ", col=rev(rainbow(10)))
邊界和網(wǎng)格
> Heatmap(heatmapdata, name=" ", col=col_fun,
+ border=TRUE,rect_gp=gpar(col="white",lwd=2))
渲染聚類樹
> #install.packages("dendextend")
> library(dendextend)
> col_dend = hclust(dist(t(heatmapdata)))
> row_dend = hclust(dist(heatmapdata))
> Heatmap(heatmapdata, name=" ",
+ cluster_columns=color_branches(col_dend,k=4),
+ cluster_rows = color_branches(row_dend, k=4))
組內(nèi)聚類
> Heatmap(heatmapdata, name=" ",
+ cluster_columns=cluster_within_group(heatmapdata,anno_col$Group))
分割熱圖
> Heatmap(heatmapdata, name=" ",
+ split = anno_row$chromosome,
+ column_split = anno_col$Group)
> Heatmap(heatmapdata, name=" ",
+ split = anno_row$chromosome,
+ column_split = anno_col$Group,
+ column_title = NULL,
+ row_title_gp = gpar(col="white",fontsize=11,
+ fill=brewer.pal(8,"Set1")[1:5]))
圖例
> Heatmap(heatmapdata, name=" ",
+ heatmap_legend_param = list(legend_height=unit(6,"cm"),
+ grid_width=unit(0.8,"cm")))
簡單注釋
> ha_col <- HeatmapAnnotation(Group=anno_col$Group,
+ col=list(Group=Group_col))
> ha_row <- HeatmapAnnotation(Chr=anno_row$chromosome,
+ col=list(Chr=Chr_col),
+ which="row")
> Heatmap(heatmapdata, name=" ",
+ top_annotation = ha_col,
+ left_annotation = ha_row)
復(fù)雜注釋
> ha_box <- HeatmapAnnotation(boxplot=anno_boxplot(heatmapdata,
+ height=unit(2,"cm"),
+ box_width = 0.9,
+ gp=gpar(fill=Chr_col[anno_row$chromosome]),
+ outline = F),
+ which="row")
> Heatmap(heatmapdata, name=" ",
+ top_annotation = ha_col,
+ left_annotation = ha_row,
+ right_annotation = ha_box)
標簽注釋
> ha_lab <- rowAnnotation(gene=anno_mark(at=which(row.names(heatmapdata)=="rna34188"),
+ labels = "CycG"))
> Heatmap(heatmapdata, name=" ",
+ show_row_names = F,
+ right_annotation = ha_lab)
(2)reference
可以去看ComplexHeatmap Complete Reference嘱能,非常詳細。