復(fù)雜圖像

我們可以使用一些簡(jiǎn)單的圖形對(duì)象圆米，將它們組合起來(lái)吊输，繪制成各種復(fù)雜的圖形

1. 圓形柱狀圖

我們介紹過(guò)了各種柱狀圖/條形圖的繪制，如堆積型激率、并列型咳燕，還有徑向柱狀圖，現(xiàn)在我們?cè)俳榻B一種圓形柱狀圖的繪制柱搜。

圓形柱狀圖就是將圓形布局中迟郎，每一條軌跡當(dāng)做一個(gè)柱子。例如聪蘸，我們有 9 個(gè)樣本宪肖，要繪制每個(gè)樣本的腫瘤純度信息

首先，構(gòu)造數(shù)據(jù)

category <- paste0("sample", "_", 1:9)
percent <- sort(sample(40:80, 9))
color <- rev(rainbow(length(percent)))

初始化健爬，以 12 點(diǎn)鐘方向作為起始控乾，同時(shí)只設(shè)置一個(gè)扇形，即整個(gè)圓就是一個(gè) sector

circos.par("start.degree" = 90, cell.padding = c(0, 0, 0, 0))
circos.initialize("a", xlim = c(0, 100))

添加圖形

circos.track(
  ylim = c(0.5, length(percent)+0.5), track.height = 0.8, 
  bg.border = NA,
  panel.fun = function(x, y) {
    xlim = CELL_META$xlim
    # 添加圓形中心線
    circos.segments(rep(xlim[1], 9), 1:9,
                    rep(xlim[2], 9), 1:9,
                    col = "#CCCCCC")
    # 添加 9 個(gè)圓形矩形
    circos.rect(rep(0, 9),
                1:9 - 0.45,
                percent,
                1:9 + 0.45,
                col = color,
                border = "white")
    # 添加文本信息
    circos.text(
      rep(xlim[1], 9),
      1:9,
      paste(category, " - ", percent, "%"),
      facing = "downward",
      adj = c(1.05, 0.5),
      cex = 0.8
    )
    # 添加軸信息
    breaks = seq(0, 85, by = 5)
    circos.axis(
      h = "top",
      major.at = breaks,
      labels = paste0(breaks, "%"),
      labels.cex = 0.6
    )
  })

2. 直方圖

可以使用 circos.trackHist() 函數(shù)在所有單元格中繪制直方圖娜遵，如果設(shè)置 draw.density = TRUE 則會(huì)繪制數(shù)據(jù)分布的密度曲線蜕衡。

默認(rèn)情況下，每個(gè)單元格的柱子的數(shù)量會(huì)根據(jù)數(shù)據(jù)自動(dòng)確定设拟，通過(guò)固定 bin.size 的值慨仿，可以讓所有的單元格繪制相同數(shù)量的柱子，有利于單元格之間的分布比較纳胧。

例如

x <- rnorm(1600)
sectors <- sample(letters[1:16], 1600, replace = TRUE)
circos.initialize(sectors, x = x)
circos.trackHist(
  sectors, x = x, col = "#a6cee3", 
  border = "#1f78b4"
  )
circos.trackHist(
  sectors, x = x, bin.size = 0.1,  
  col = "#b2df8a", border = "#33a02c"
  )
circos.trackHist(
  sectors, x = x, draw.density = TRUE, 
  col = "#fdbf6f", border = "#ff7f00"
  )
circos.clear()

3. 系統(tǒng)發(fā)育樹(shù)

圓形樹(shù)狀圖可以有多方面的應(yīng)用镰吆，比如系統(tǒng)發(fā)育樹(shù)的展示。

R 中有多種樹(shù)結(jié)構(gòu)類(lèi)跑慕，例如 hclust, dendrogram 和 phylo万皿，它們之間可以進(jìn)行相互轉(zhuǎn)換摧找，所以我們只使用 dendrogram 類(lèi)來(lái)說(shuō)明

首先，導(dǎo)入 ape 包提供的鳥(niǎo)類(lèi)物種信息牢硅，并對(duì)數(shù)據(jù)進(jìn)行層次聚類(lèi)

library(ape)

data(bird.orders)
hc <- as.hclust(bird.orders)

然后蹬耘，使用 cutree() 將數(shù)進(jìn)行切割，劃分出 6 個(gè)物種

# 獲取鳥(niǎo)類(lèi)名稱
labels <- hc$labels
# 劃分為 6 個(gè)物種
ct <- cutree(hc, 6)
# 物種的數(shù)量
n <- length(labels) 
# 轉(zhuǎn)換為 dendrogram 結(jié)構(gòu)
dend <- as.dendrogram(hc)

因?yàn)闃?shù)形圖在內(nèi)側(cè)减余，所以要先繪制物種的標(biāo)簽

circos.par(cell.padding = c(0, 0, 0, 0))
# 只需要一個(gè)扇形
circos.initialize("a", xlim = c(0, n))
circos.track(
  ylim = c(0, 1), bg.border = NA, track.height = 0.3, 
  panel.fun = function(x, y) {
    for (i in seq_len(n)) {
      circos.text(
        i - 0.5, 0, labels[i], adj = c(0, 0.5),
        facing = "clockwise", niceFacing = TRUE,
        col = ct[labels[i]], cex = 0.5
      )
    }
  })

最后使用 circos.dendrogram() 函數(shù)繪制樹(shù)形圖

library(dendextend)

dend <- color_branches(dend, k = 6, col = 1:6)
dend_height <- attr(dend, "height")
circos.track(
  ylim = c(0, dend_height), bg.border = NA, 
  track.height = 0.4, 
  panel.fun = function(x, y) {
    circos.dendrogram(dend)
    }
  )
circos.clear()

可以使用 dendextend 包對(duì)樹(shù)形圖的屬性進(jìn)行設(shè)置综苔，樹(shù)形圖默認(rèn)是朝外的，可以設(shè)置 facing = "inside"佳励，使其朝向內(nèi)部

circos.dendrogram(dend, facing = "inside")

注意休里，要先繪制樹(shù)形圖后繪制標(biāo)簽

4. 手動(dòng)繪制圓形熱圖

既然我們可以繪制圓形樹(shù)狀圖，那么很容易想到赃承，在最外層使用 circos.rect() 函數(shù)添加一圈熱圖妙黍，就變成了圓形熱圖。

我們要繪制兩個(gè)獨(dú)立的熱圖瞧剖，首先拭嫁，構(gòu)造數(shù)據(jù)

mat <- matrix(rnorm(100*10), nrow = 100, ncol = 10)
col_fun <- colorRamp2(c(-2, 0, 2), c("#fc8d59", "#ffffbf", "#99d594"))
# 設(shè)置兩個(gè)扇形
sectors <- rep(letters[1:2], times = c(30, 70))
mat_list <- list(
  a = mat[sectors == "a", ],
  b = mat[sectors == "b", ]
  )
# 將聚類(lèi)結(jié)果轉(zhuǎn)換為 dendrogram 類(lèi)
dend_list <- list(
  a = as.dendrogram(hclust(dist(mat_list[["a"]]))),
  b = as.dendrogram(hclust(dist(mat_list[["b"]])))
)

我們需要從外到內(nèi)依次繪制，先繪制圓形熱圖

circos.par(cell.padding = c(0, 0, 0, 0), gap.degree = 5)
circos.initialize(sectors, xlim = cbind(c(0, 0), table(sectors)))
circos.track(
  ylim = c(0, 10), bg.border = NA, 
  panel.fun = function(x, y) {
    sector.index = CELL_META$sector.index
    m = mat_list[[sector.index]]
    dend = dend_list[[sector.index]]
    
    m2 = m[order.dendrogram(dend),]
    col_mat = col_fun(m2)
    nr = nrow(m2)
    nc = ncol(m2)
    for (i in 1:nc) {
      circos.rect(
        1:nr - 1, rep(nc - i, nr), 1:nr,
        rep(nc - i + 1, nr), border = col_mat[, i],
        col = col_mat[, i])
    }
})

因?yàn)槲覀冃枰L制了兩個(gè)樹(shù)狀圖抓于，要保證它們的高度一致做粤，可以去兩個(gè)當(dāng)中更高的那個(gè)來(lái)設(shè)置 ylim

# 獲取最大樹(shù)高
max_height <- max(sapply(dend_list, function(x) attr(x, "height")))
circos.track(
  ylim = c(0, max_height), bg.border = NA, 
  track.height = 0.3, panel.fun = function(x, y) {
    sector.index = get.cell.meta.data("sector.index")
    dend = dend_list[[sector.index]]
    circos.dendrogram(dend, max_height = max_height)
  }
)
circos.clear()

圖例

circlize 為用戶提供了完全自由的圖形設(shè)計(jì)，但是缺少了對(duì)圖例的控制捉撮。

例如怕品，對(duì)于如下圖形

col_fun <- colorRamp2(c(-2, 0, 2), c("green", "blue", "red"))
circlize_plot <- function() {
  set.seed(123)
  sectors = letters[1:10]
  circos.initialize(sectors, xlim = c(0, 1))
  circos.track(ylim = c(0, 1), panel.fun = function(x, y) {
    circos.points(runif(20), runif(20), cex = 0.5, pch = 16, col = 2)
    circos.points(runif(20), runif(20), cex = 0.5, pch = 16, col = 3)
  })
  circos.track(ylim = c(0, 1), panel.fun = function(x, y) {
    circos.lines(sort(runif(20)), runif(20), col = 3)
    circos.lines(sort(runif(20)), runif(20), col = 7)
  })
  
  for(i in 1:10) {
    circos.link(sample(sectors, 1), sort(runif(10))[1:2], 
                sample(sectors, 1), sort(runif(10))[1:2],
                col = add_transparency(col_fun(rnorm(1))))
  }
  circos.clear()
}
circlize_plot()

現(xiàn)在要為這三個(gè)軌跡添加圖例，我們可以使用 ComplexHeatmap 的 Legend() 函數(shù)來(lái)自定義圖例

library(ComplexHeatmap)
# 點(diǎn)
lgd_points <- Legend(
  at = c("label1", "label2"), type = "points", 
  legend_gp = gpar(col = 2:3), title_position = "topleft", 
  title = "Track1"
)
# 線
lgd_lines <- Legend(
  at = c("label3", "label4"), type = "lines", 
  legend_gp = gpar(col = 4:5, lwd = 2), 
  title_position = "topleft", title = "Track2"
)
# 顏色條
lgd_links <- Legend(
  at = c(-2, -1, 0, 1, 2), col_fun = col_fun, 
  title_position = "topleft", title = "Links"
)

要將這三個(gè)圖例合并在一起巾遭，可以使用 packLegend() 函數(shù)肉康，默認(rèn)按豎直方向添加

lgd_list_vertical <- packLegend(lgd_points, lgd_lines, lgd_links)

可以使用 draw() 函數(shù)來(lái)繪制圖例

draw(
  lgd_list_vertical, x = unit(5, "mm"), 
  y = unit(4, "mm"), just = c("left", "bottom")
)

雖然 circlize 是基于基礎(chǔ)圖形系統(tǒng)，而 ComplexHeatmap 是基于 grid 繪圖系統(tǒng)灼舍，但是兩種系統(tǒng)可以混合使用吼和。事實(shí)上，它們是繪制在同一個(gè)圖形設(shè)備上的不同圖層

但是這種方法很容易出現(xiàn)圖例與圓形布局重疊的情況骑素，更好的方式是將圖例和圓形圖分為兩部分炫乓。

我們可以用 grid 包來(lái)對(duì)繪圖區(qū)域進(jìn)行重排，例如

grid.newpage()
circle_size = unit(1, "snpc")
pushViewport(
  viewport(
    x = 0, y = 0.5, width = circle_size, 
    height = circle_size, just = c("left", "center")
    )
  )
par(omi = c(0, 0, 0.5, 0.5), new = TRUE)
circlize_plot()
upViewport()
draw(lgd_list_vertical, x = circle_size, just = "left")

使用 par(new = TRUE) 防止基礎(chǔ)圖形繪制到新的圖片上献丑，并用 omi 參數(shù)來(lái)設(shè)置邊距末捣。

也可以設(shè)置水平排列的圖例

lgd_points <- Legend(
  at = c("label1", "label2"), type = "points", 
  legend_gp = gpar(col = 2:3), title_position = "topleft",
  title = "Track1", nrow = 1
)

lgd_lines <- Legend(
  at = c("label3", "label4"), type = "lines", 
  legend_gp = gpar(col = 4:5, lwd = 2), title_position = "topleft",
  title = "Track2", nrow = 1
)

lgd_links <- Legend(
  at = c(-2, -1, 0, 1, 2), col_fun = col_fun, 
  title_position = "topleft", title = "Links", direction = "horizontal"
  )

lgd_list_horizontal <- packLegend(
  lgd_points, lgd_lines, lgd_links, 
  direction = "horizontal"
)

plot.new()
pushViewport(viewport(
  x = 0.5, y = 1, width = circle_size, 
  height = circle_size, just = c("center", "top"))
  )
par(omi = c(0, 0, 0, 0), new = TRUE)
circlize_plot()
upViewport()

draw(lgd_list_horizontal, y = unit(1, "npc") - circle_size, just = "top")

這種方式也是要調(diào)整繪圖區(qū)域的大小，來(lái)顯示圖例