隨著生物學背景知識的增加，單細胞圖譜的可視化直接用10X的Loup或者seurat的Dimplot函數(shù)直接繪制的umap/tsne圖往往很難達到要求了，這就要求我們提高繪圖技能。我們都知道ggplot2是一款很好的繪圖R包，甚至可以說在語法上擴展了R語言本身。那么，當我們需要繪圖的時候，自然我們會想到它及其周邊。今天我們就主要地看一下ggforce這個包帶給我們的可能性。

首先，我載入數(shù)據(jù)：

library(Seurat)
library(ggplot2)
library(tidyverse)
pbmc <- readRDS('G:\\Desktop\\Desktop\\RStudio\\single_cell\\filtered_gene_bc_matrices\\hg19pbmc_tutorial.rds')
pbmc


An object of class Seurat 
13714 features across 2638 samples within 1 assay 
Active assay: RNA (13714 features)
 3 dimensional reductions calculated: pca, umap, tsne···

為了提供更多的分群結果，我們再跑一次FindClusters,當然你也可自己構建分組方式，比如不同的樣本，VDJ不同的克隆型等

pbmc<-FindClusters(pbmc,resolution = c(.4,.8,1.2,1.6,2)) 
head(pbmc@meta.data)

               orig.ident nCount_RNA nFeature_RNA percent.mt RNA_snn_res.0.5
AAACATACAACCAC     pbmc3k       2419          779  3.0177759               1
AAACATTGAGCTAC     pbmc3k       4903         1352  3.7935958               3
AAACATTGATCAGC     pbmc3k       3147         1129  0.8897363               1
AAACCGTGCTTCCG     pbmc3k       2639          960  1.7430845               2
AAACCGTGTATGCG     pbmc3k        980          521  1.2244898               6
AAACGCACTGGTAC     pbmc3k       2163          781  1.6643551               1
               seurat_clusters RNA_snn_res.0.4 RNA_snn_res.0.8 RNA_snn_res.1.2
AAACATACAACCAC               1               2               1               5
AAACATTGAGCTAC               0               3               2               2
AAACATTGATCAGC               2               2               1               1
AAACCGTGCTTCCG               3               1               4               4
AAACCGTGTATGCG               8               6               7               8
AAACGCACTGGTAC               1               2               1               1
               RNA_snn_res.1.6 RNA_snn_res.2
AAACATACAACCAC               9             1
AAACATTGAGCTAC               2             0
AAACATTGATCAGC               1             2
AAACCGTGCTTCCG               4             3
AAACCGTGTATGCG               8             8
AAACGCACTGGTAC               1             1

為了調用ggplot2我們把UMAP的坐標放到metadata中：

pbmc<-AddMetaData(pbmc,pbmc@reductions$umap@cell.embeddings,col.name = colnames(pbmc@reductions$umap@cell.embeddings))

head(pbmc@meta.data)

讀入一套我珍藏多年的顏色列表：

allcolour=c("#DC143C","#0000FF","#20B2AA","#FFA500","#9370DB","#98FB98","#F08080","#1E90FF","#7CFC00","#FFFF00",
              "#808000","#FF00FF","#FA8072","#7B68EE","#9400D3","#800080","#A0522D","#D2B48C","#D2691E","#87CEEB","#40E0D0","#5F9EA0",
              "#FF1493","#0000CD","#008B8B","#FFE4B5","#8A2BE2","#228B22","#E9967A","#4682B4","#32CD32","#F0E68C","#FFFFE0","#EE82EE",
              "#FF6347","#6A5ACD","#9932CC","#8B008B","#8B4513","#DEB887")

我用ggplot畫一個帶有標簽的umap圖：

class_avg <- pbmc@meta.data %>%
    group_by(RNA_snn_res.2) %>%
    summarise(
      UMAP_1 = median(UMAP_1),
      UMAP_2 = median(UMAP_2)
    )

    
umap <-  ggplot(pbmc@meta.data ,aes(x=UMAP_1,y=UMAP_2))+
    geom_point(aes(color=RNA_snn_res.2))+
    scale_color_manual(values = allcolour)+
    geom_text(aes(label = RNA_snn_res.2), data = class_avg)+
    theme(text=element_text(family="Arial",size=18)) +
    theme(panel.background = element_rect(fill='white', colour='black'), 
                     panel.grid=element_blank(), axis.title = element_text(color='black',
                                                      family="Arial",size=18),axis.ticks.length = unit(0.4,"lines"), 
                     axis.ticks = element_line(color='black'), 
                     axis.ticks.margin = unit(0.6,"lines"),
                     axis.line = element_line(colour = "black"), 
                     axis.title.x=element_text(colour='black', size=18),
                     axis.title.y=element_text(colour='black', size=18),
                     axis.text=element_text(colour='black',size=18),
                     legend.title=element_blank(),
                     legend.text=element_text(family="Arial", size=18),
                     legend.key=element_blank())+
    theme(plot.title = element_text(size=22,colour = "black",face = "bold"))  + 
    guides(colour = guide_legend(override.aes = list(size=5)))

umap

學會了嗎?

為了使我們的圖層不要那么復雜，還是先畫一個簡單的：

umap <-  ggplot(pbmc@meta.data ,aes(x=UMAP_1,y=UMAP_2,color=RNA_snn_res.2))+
  geom_point() 

umap

好奇的我們來看一下umap這個圖都有什么：

head(umap$data)
umap$theme
umap$layers
umap$scales
umap$mapping
umap$coordinates
umap$facet
umap$plot_env  
umap$labels

然后，我們請出ggforce這個包，看看第一次的驚喜。

library(ggforce)  
umap + facet_zoom(x = RNA_snn_res.2 == "14")

想要細致刻畫某個亞群，這不失是一個方法：

umap +  facet_zoom(xlim = c(-15, -10), ylim = c(0, 2.5))

umap + geom_mark_rect(aes(label = RNA_snn_res.2), show.legend = FALSE) +
  theme_void() 

給每個群加框加標簽，優(yōu)雅：

library(concaveman)
umap +
geom_mark_hull(aes(label = RNA_snn_res.2)) +
  theme_void()  

可以根據(jù)自己的數(shù)據(jù)格式換一換 ：

umap +  
geom_mark_hull(aes(label = RNA_snn_res.2, fill = RNA_snn_res.0.4), show.legend = FALSE) +
  theme_void() 

如果有需要特別化為一類的可以用背景色來圈?。?/p>

umap +  
geom_mark_hull(aes(label = RNA_snn_res.2, fill = RNA_snn_res.2), show.legend = FALSE, expand = unit(3, "mm")) +
  theme_void() 

umap +  
  geom_mark_hull(aes(label = RNA_snn_res.2, fill = RNA_snn_res.2), show.legend = FALSE, expand = unit(3, "mm")) +
  theme_no_axes() 

desc <- 'I am a  luck dog'
umap +
  geom_mark_ellipse(aes(filter = RNA_snn_res.2 == '14', label = '14', 
                        description = desc))  

想對某一亞群做進一步的注釋：

你也可以：

umap +
  geom_mark_hull(aes(filter = RNA_snn_res.2 == '14', label = '14', 
                        description = desc)) 

umap +
geom_voronoi_tile(aes(fill = RNA_snn_res.2, group = -1L)) + 
  geom_voronoi_segment()

Voronoi圖背后的想法是將圖的區(qū)域分割成盡可能多的部分。與網(wǎng)格或熱圖不同，Voronoi根據(jù)與其他點的接近程度為每個點繪制自定義形狀。它返回一個看起來像彩色玻璃的圖。這可以很好地確定每個區(qū)域內的最近點。例如，零售商可以使用它來查看他們的商店位置所覆蓋的區(qū)域，并可以幫助他們做出決策，根據(jù)每個Voronoi形狀的大小來優(yōu)化他們的位置。

其實這并看不清

umap +
geom_voronoi_tile(aes(fill = RNA_snn_res.2), max.radius = 1,colour = 'black')

看一看出哪些地方的細胞比較密集，這一點當然需要好的降維結構了，細胞密集與否分別代表什么？越密集的區(qū)域細胞距離越近，說明異質性較低。當然，這和降維結構有關。

umap +
geom_voronoi_tile(aes(fill = RNA_snn_res.2), max.radius = .1,colour = 'black')

umap + 
  geom_mark_hull(aes(fill = RNA_snn_res.2), expand = unit(3, "mm")) +
  coord_cartesian(xlim = c(-15, -10), ylim = c(0, 2.5))  +
  geom_voronoi_segment()

有種細胞的感覺了嗎?

最后，作為附贈：

pbmc@meta.data %>%
  gather_set_data(6:11) %>%
  ggplot(aes(x, id = id, split = y, value = 1))  +
  geom_parallel_sets(aes(fill = RNA_snn_res.0.4), show.legend = FALSE, alpha = 0.3) +
  geom_parallel_sets_axes(axis.width = 0.1, color = "lightgrey", fill = "white") +
  geom_parallel_sets_labels(angle = 0) +
  theme_no_axes()

pbmc@meta.data %>%
  count(RNA_snn_res.0.4) %>%
  ggplot() +
  geom_arc_bar(aes(x0 = 0, y0 = 0, r0 = 0.7, r = 1, amount = n, fill = RNA_snn_res.0.4), alpha = 0.3, stat = "pie") 

p1 <- pbmc@meta.data %>%
  count(RNA_snn_res.0.4) %>%
  mutate(focus = ifelse(RNA_snn_res.0.4 == "0", 0.2, 0)) %>%
  ggplot()+
  geom_arc_bar(aes(x0 = 0, y0 = 0, r0 = 0.7, r = 1, amount = n, 
                   fill = RNA_snn_res.0.4, explode = focus), 
               alpha = 1, stat = "pie") +
  scale_fill_manual(values = allcolour)

• sc-RAN-seq 數(shù)據(jù)分析||Seurat新版教程:Guided Clustering Tutorial
• sc-RAN-seq 數(shù)據(jù)分析||Seurat新版教程: Integrating datasets to learn cell-type specific responses
• sc-RAN-seq 數(shù)據(jù)分析||Seurat新版教程: Using sctransform in Seurat
• 單細胞轉錄組數(shù)據(jù)分析||Seurat新版教程:Differential expression testing
• 單細胞轉錄組 數(shù)據(jù)分析||Seurat新版教程:New data visualization methods in v3.0
• 單細胞轉錄組數(shù)據(jù)分析||Seurat并行策略
• Seurat Weekly NO.0 || 開刊詞
• Seurat Weekly NO.1 || 到底分多少個群是合適的？！
• Seurat Weekly NO.2 || 我該如何取子集
• 你到底想要什么樣的umap/tsne圖？
• scRNA-seq擬時分析 || Monocle2 踩坑教程
• scRNA-seq數(shù)據(jù)分析 || Monocle3

https://www.r-bloggers.com/the-ggforce-awakens-again/
https://rviews.rstudio.com/2019/09/19/intro-to-ggforce/