前言

Hello小伙伴們大家好，我是生信技能樹的小學徒”我才不吃蛋黃“。今天是肺腺癌單細胞數據集GSE189357復現系列第三期。第二期我們對細胞亞群進行了注釋（肺腺癌單細胞數據集GSE189357復現(二)：細胞注釋）。

本期，我們將使用多種可視化方法，繪制FeaturePlot，ggplot，DoHeatmap圖。

1.背景介紹

在胃癌系列推文中，我提到繪好看圖需要審美能力強。實際上，繪圖能力也能很好的反映我們的科研素養(yǎng)。繪圖能力與科研能力之間的關系非常緊密，因為有效的可視化在科學研究中扮演著至關重要的角色。繪圖不僅是展示研究成果的手段，更是理解復雜數據、發(fā)現新知識以及向同行清晰傳達信息的核心工具。在單細胞數據可視化中，我們大部分情況下都是充當“調包俠”，使用各種R包繪制各種圖形。

R語言是一種廣泛使用的統(tǒng)計分析和圖形表示語言，它在數據可視化方面有著強大的功能，但是也有其缺點。

R語言可視化的優(yōu)點有很多：

豐富的圖形庫：R擁有大量的圖形庫，如ggplot2、lattice、base plots等，這些庫提供了豐富的圖形類型和定制選項。

高度可定制：R的圖形可以高度定制，從顏色、形狀到布局，用戶可以精細控制圖形的每一個方面。

數據操作能力強：R是一種強大的數據分析工具，它可以在繪圖之前對數據進行復雜的處理和分析。

集成統(tǒng)計分析：R語言在統(tǒng)計分析方面非常強大，可以直接在繪圖中集成統(tǒng)計測試和模型結果。

動態(tài)報告：R的可視化可以輕松地集成到動態(tài)報告中，如R Markdown和Shiny，這使得結果可以交互式地呈現。

社區(qū)支持：R有一個活躍的社區(qū)，用戶可以從中獲得大量的教程、論壇討論和包更新。

開源免費：R是開源軟件，可以免費使用和修改，這使得它在學術界和數據科學領域非常受歡迎。

跨平臺：R可以在多種操作系統(tǒng)上運行，包括Windows、macOS和Linux。

R語言缺點：

學習曲線：對于初學者來說，R的學習曲線可能比較陡峭，尤其是對于那些沒有編程背景的用戶。

性能問題：雖然R在數據處理和可視化方面表現出色，但它在處理非常大的數據集時可能會變慢。

圖形更新：在R中更新圖形可能需要重新運行整個腳本，這在迭代過程中可能會顯得繁瑣。

3D可視化有限：雖然R可以創(chuàng)建3D圖形，但與一些專門的軟件相比，它的3D可視化能力相對有限。

用戶界面：R的默認用戶界面（RStudio除外）可能不如一些集成開發(fā)環(huán)境（IDE）友好。

依賴管理：R的包依賴管理有時可能會導致版本沖突，需要用戶手動解決。

專業(yè)圖形制作：雖然R可以制作高質量的圖形，但與專業(yè)的圖形設計軟件相比，它在制作出版級別的圖形方面可能不夠直觀。

交互性：雖然R可以通過Shiny等工具創(chuàng)建交互式圖形，但這些交互性可能不如專門的交互式可視化工具（如Tableau）直觀和強大。

R繪圖最重要的一點是無法像PS一樣隨意改變圖形，只能在參數限定范圍內修改調整。因此，要想把圖做的好看，就要選擇合適的繪圖函數，并不停的調參。所以說，理解函數是R繪圖的基礎。

2.可視化

首先加載R包，創(chuàng)建新的文件夾，讀取細胞注釋后的Seurat數據：

rm(list=ls())
library(Seurat)
library(ggplot2)
library(clustree)
library(cowplot)
library(dplyr)
library(SingleR)
library(celldex)
#BiocManager::install("celldex")
# library(devtools)
# install_github("arc85/singleseqgset")
library(singleseqgset)
library(devtools)
library(grid)
library(gridExtra)
getwd()
dir.create("4-plot")
setwd('4-plot/')
sce.all=readRDS( "../3-Celltype/sce_celltype.rds")
sce.all
#Idents(sce.all)

2.1 Featureplot可視化基因

#EPCAM,NKG7,LYZ,CD79A,CLDN5,DCN
marker <- c('EPCAM','NKG7','LYZ','CD79A','CLDN5','DCN')
gene = marker
FeaturePlot(sce.all,features = marker,cols = c("lightgrey" ,"#DE1F1F"),ncol=3,raster=FALSE)
ggsave('FeaturePlot_marker.pdf',width = 12,height = 8)

我們可以看到，上圖各個基因的顏色閾值范圍不一致，我們可以調參，使其保持一致：

p1 <- FeaturePlot(sce.all, features = marker, combine = FALSE,ncol=3,raster=FALSE )
#colours = c('lightgrey', "#DE1F1F")
fix.sc <- scale_color_gradientn( colours = c('lightgrey', "#DE1F1F"),  limits = c(0, 6))
#+NoLegend()+NoAxes()
p2 <- lapply(p1, function (x) x + fix.sc)
CombinePlots(p2)

批量畫基因,注意圖例的范圍不同：

FeaturePlot(sce.all, features =marker, 
            cols = c("lightgrey", 'red'),
            ncol = 3 ) & NoLegend() & NoAxes() & theme(
              panel.border = element_rect(color = "black", size = 1)
            )

Featureplot還可以把兩個基因畫在同一個圖中,看右上角可以發(fā)現黃色越深的地方兩個基因疊加越多：

FeaturePlot(sce.all, features = c('S100A9','S100A8'),
            cols = c("lightgrey", "green", "orange"),
            blend=T,blend.threshold=0)

結合ggplot函數，我們還可以把三個基因畫在同一個圖中：

提取tsne坐標，并提取基因表達數據并與tsne坐標合并：

tsne_df <- as.data.frame(sce.all@reductions$umap@cell.embeddings)
tsne_df$cluster <- as.factor(sce.all$celltype)
head(tsne_df)
gene_df <- as.data.frame(GetAssayData(object = sce.all, slot = "data")[c('S100A9','S100A8','CXCL8'), ])

ggplot繪制圖形：

library(ggnewscale)
merged_df <- merge(t(gene_df), tsne_df, by = 0, all = TRUE)
head(merged_df)
colnames(merged_df)
ggplot(merged_df, vars = c("umap_1", "umap_2", 'S100A9','S100A8','CXCL8'), aes(x = umap_1, y = umap_2, colour = S100A9)) +
  geom_point(size=0.3, alpha=1) +
  scale_colour_gradientn(colours = c("lightgrey", "green"), limits = c(0, 0.3), oob = scales::squish) +
  new_scale_color() +
  geom_point(aes(colour = S100A8), size=0.3, alpha=0.7) +
  scale_colour_gradientn(colours = c("lightgrey", "blue"), limits = c(0.1, 0.2), oob = scales::squish) +
  new_scale_color() +
  geom_point(aes(colour = CXCL8), size=0.3, alpha=0.1) +
  scale_colour_gradientn(colours = c("lightgrey", "red"), limits = c(0, 0.3), oob = scales::squish)+
  theme_classic()

上面的圖花里胡哨，目的是什么？個人認為，可視化的主要目的是清晰的呈現數據的特征并揭示生物學意義。使用FeaturePlot將多個基因的表達繪制在同一個圖中，我們可以比較這不同基因在細胞中的共表達模式，以及探索它們在特定細胞類型或亞群中的潛在生物學相關性。

2.2 DoHeatmap繪制熱圖

DoHeatmap 在單細胞 RNA 測序數據分析中用于可視化細胞群體或類型的基因表達模式，其目的包括展示差異表達基因、識別細胞亞群、展示標記基因、比較基因表達模式、分析細胞類型的群體特征以及探索潛在的功能和通路。這種可視化工具是解釋單細胞數據、理解細胞異質性以及揭示潛在生物學機制的強大手段。

利用DoHeatmap函數繪制熱圖，可以展示不同細胞類型的top5 maker：

Idents(sce.all)
table(sce.all$celltype)
Idents(sce.all) = sce.all$celltype
sce1 = sce.all[, sce.all$celltype %in% c( 'B', 'Endothelial','Epithelial', 'Fibro','Myeloid' ,'T&NK' )]

if (!file.exists('sce.markers.csv')) {
  sce.markers <- FindAllMarkers(object = sce1, only.pos = TRUE, 
                                min.pct = 0.25, 
                                thresh.use = 0.25)
  write.csv(sce.markers,file='sce.markers.csv')
} else {
  
  sce.markers = read.csv('sce.markers.csv',row.names = 1)
}


library(dplyr) 
top5 <- sce.markers%>% group_by(cluster) %>% top_n(5, avg_log2FC)

為了防止數據量太大不好出圖，這里在每個亞群提取出來100個：

sce.Scale <- ScaleData(subset(sce1,downsample=100),
                       features = top5$gene )
DoHeatmap(sce.Scale,
          features = top5$gene ,
          # group.by = "celltype",
          assay = 'RNA', label = T)+
  scale_fill_gradientn(colors = c("white","grey","firebrick3"))

ggsave('markers_heatmap.pdf',width = 10,height = 7)

2.3 DotPlot繪制氣泡圖

top5_dotplot <- DotPlot(sce.all, features = top5$gene)+
  theme(axis.text.x = element_text(angle = 90, vjust = 1, hjust = 1))
top5_dotplot
ggsave('markers_top5_dotplot.pdf',width = 10,height = 7)

setwd('../')

3.參數解析

本期單細胞基因可視化用到的三種函數分別為FeaturePlot()、DoHeatmap()和DotPlot()。代碼參數及解析如下：

FeaturePlot()：

FeaturePlot(
object,   #Seurat對象，包含要進行可視化的數據集。
features,  #字符向量，指定要在圖中展示的特征（例如基因或元數據列名）。
dims = c(1, 2),  #數字向量，長度為2，指定要用于繪圖的降維空間的維度（例如c(1, 2)表示第一和第二維度）。
cells = NULL,  #可選的，指定要在圖中展示的細胞的向量。默認是所有細胞。
cols = if (blend) { c("lightgrey", "#ff0000", "#00ff00") } else {
c("lightgrey", "blue") },
pt.size = NULL,
order = FALSE,
min.cutoff = NA,
max.cutoff = NA,
reduction = NULL,
split.by = NULL,
shape.by = NULL,
slot = "data",
blend = FALSE,
blend.threshold = 0.5,
label = FALSE,
label.size = 4,
repel = FALSE,
ncol = NULL,
coord.fixed = FALSE,
by.col = TRUE,
sort.cell = NULL,
interactive = FALSE,
combine = TRUE
)

FeaturePlot函數參數解析：

object: Seurat對象，包含要進行可視化的數據集。

features: 字符向量，指定要在圖中展示的特征（例如基因或元數據列名）。

dims: 數字向量，長度為2，指定要用于繪圖的降維空間的維度（例如c(1, 2)表示第一和第二維度）。

cells: 可選的，指定要在圖中展示的細胞的向量。默認是所有細胞。

cols: 字符向量或顏色向量，用于定義繪圖中使用的顏色漸變。

pt.size: 點的大小。

order: 邏輯值，指定是否根據特征表達量的順序繪制細胞，可以幫助突出表達某特征的細胞。

min.cutoff, max.cutoff: 用于指定每個特征的表達值截止的向量?？梢允褂冒俜治粩祦碇付ń刂裹c（例如，'q1', 'q99'表示1%和99%分位數）。

reduction: 指定要使用的降維技術，例如"umap"、"tsne"或"pca"。如果未指定，FeaturePlot會依次查找"umap"、"tsne"、"pca"中可用的結果。

split.by: 字符串，指定元數據中的一個變量，用于按該變量的不同類別拆分并分別繪制圖形。

shape.by: 可選的，允許根據某個細胞屬性改變點的形狀。

slot: 指定從哪個Seurat對象的槽中提取表達數據進行可視化。

blend: 邏輯值，指定是否混合兩個特征的表達值來同時可視化它們。

blend.threshold: 設置用于混合特征表達的閾值，范圍從0到1。

label: 是否在圖上標記群組。

label.size: 標簽文字的大小。

repel: 邏輯值，指定是否使用標簽排斥機制，以避免標簽之間的重疊。

ncol: 數字，指定在使用split.by時，合并到一個圖中的列數。

coord.fixed: 邏輯值，指定是否使用固定的縱橫比繪制坐標系。

by.col: 邏輯值，指定在分列展示時是否按列而非按行展示特征圖。

interactive: 邏輯值，指定是否啟用交互式FeaturePlot。

combine: 邏輯值，指定是否將多個特征圖合并為單個繪圖對象。如果為FALSE，則返回一個包含多個ggplot對象的列表。

DoHeatmap()

DoHeatmap(
object,
features = NULL,
cells = NULL,
group.by = "ident",
group.bar = TRUE,
group.colors = NULL,
disp.min = -2.5,
disp.max = NULL,
slot = "scale.data",
assay = NULL,
label = TRUE,
size = 5.5,
hjust = 0,
angle = 45,
raster = TRUE,
draw.lines = TRUE,
lines.width = NULL,
group.bar.height = 0.02,
combine = TRUE
)

DoHeatmap()函數參數解析：

object :  一個Seurat對象，包含要進行可視化的數據集。

features : 要打印的特征向量，默認為variableffeatures（object=object）

cells : 要繪制的細胞向量

group.by : 一種變量向量，用于按單元格分組；將“ident”傳遞給按單元格分組的標識類

group.bar : 添加顯示單元格組狀態(tài)的顏色欄

group.colors : 要用于顏色欄的顏色

disp.min : 最小顯示值（以下所有值均被剪裁）

disp.max : 最大顯示值（以上所有值均被剪裁）；如果插槽為'比例數據，否則為6

slot : 要使用的數據槽，從'原始數據'、'數據'或'比例數據'

assay : 從中提取

label : 在顏色欄上方標記單元格標識

size : 顏色欄上方文本的大小

hjust : 顏色欄上方文本的水平對正

angle : 顏色欄上方文本的角度

raster : 如果為真，則使用//oraprdnt繪制，否則使用available。//oraprdnt在某些查看應用程序（如預覽）上，由于光柵是如何插值的，因此可能會顯得模糊。如果遇到該問題，請將其設置為false（請注意，打印可能需要更長的時間來生成/渲染）。

draw.lines : 包括分隔組的白線

lines.width : 整數，用于調整分隔白色的寬度行。對應每個組之間的“細胞”數量。

group.bar.height : 縮放顏色欄的高度

combine : 將繪圖合并到單個PatchworkedgPlot對象中。如果為FALSE，則返回ggplot對象的列表

DotPlot()

DotPlot(
object,
assay = NULL,
features,
cols = c("lightgrey", "blue"),
col.min = -2.5,
col.max = 2.5,
dot.min = 0,
dot.scale = 6,
idents = NULL,
group.by = NULL,
split.by = NULL,
cluster.idents = FALSE,
scale = TRUE,
scale.by = "radius",
scale.min = NA,
scale.max = NA
)

DotPlot()函數參數解析：

object : 一個Seurat對象，包含要進行可視化的數據集。

assay : 要使用的分析名稱，默認為活動分析

features : 特征的輸入向量，或特征向量的命名列表如果需要特征分組面板（復制舊SplitDotPlotGG的功能）

cols : 要打印的顏色：來自rcolorbrewer:的調色板的名稱：布魯爾.pal.info，一對定義漸變的顏色，或3+個定義多個漸變的顏色（如果拆分方式（已設置）

col.min : 最小縮放平均表達式閾值（everythingsmaller將設置為該值）

col.max : 最大縮放平均表達式閾值（所有較大的值都將設置為該值）

dot.min : 繪制最小點的單元格分數（默認值為0）。所有表達givengene的細胞組都沒有畫點。

dot.scale : 縮放點的大小，類似于cex

idents : 要包含在繪圖中的標識類（默認為all）

group.by : 將細胞分組的因子

split.by : 用于拆分組的因子（復制舊的SplitDotPlotGG的功能）；有關詳細信息，請參閱FetchData

cluster.idents : 是否根據給定的特征按層次聚類對標識進行排序，默認值為FALSE

scale : 確定數據是否縮放，默認為TRUE

scale.by : 按“大小”或“半徑”縮放點的大小

scale.min : 設置縮放下限，默認使用NA

scale.max : 設置縮放上限，默認使用NA

結語

本期，我們使用多種可視化方法，繪制DimPlot，FeaturePlot，ggplot，DoHeatmap圖，展示了不同細胞類型基因的表達。單細胞數據分析可視化的方法遠不止以上幾種，新的繪圖R包層出不窮，各種好看的圖形讓我們目不暇接。在不停的學習代碼之余，我們還要清晰的認識到，可視化的主要目的是清晰的呈現數據的特征并揭示生物學意義。所以，我們還是要不停的看教材、看文獻、學統(tǒng)計，不停的思考。任重而道遠，我們仍需努力。下一期，我們將在此基礎上，繪制餅圖、堆積柱狀圖、箱線圖、氣泡圖等，比較不同分組之間細胞比例差異。干貨滿滿，歡迎大家持續(xù)追更，謝謝！