肺腺癌单细胞数据集GSE189357复现(三)：基因可视化

前言

Hello小伙伴们大家好，我是生信技能树的小学徒”我才不吃蛋黄“。今天是肺腺癌单细胞数据集GSE189357复现系列第三期。第二期我们对细胞亚群进行了注释（肺腺癌单细胞数据集GSE189357复现(二)：细胞注释）。

本期，我们将使用多种可视化方法，绘制FeaturePlot，ggplot，DoHeatmap图。

1.背景介绍

在胃癌系列推文中，我提到绘好看图需要审美能力强。实际上，绘图能力也能很好的反映我们的科研素养。绘图能力与科研能力之间的关系非常紧密，因为有效的可视化在科学研究中扮演着至关重要的角色。绘图不仅是展示研究成果的手段，更是理解复杂数据、发现新知识以及向同行清晰传达信息的核心工具。在单细胞数据可视化中，我们大部分情况下都是充当“调包侠”，使用各种R包绘制各种图形。

R语言是一种广泛使用的统计分析和图形表示语言，它在数据可视化方面有着强大的功能，但是也有其缺点。

R语言可视化的优点有很多：

丰富的图形库：R拥有大量的图形库，如ggplot2、lattice、base plots等，这些库提供了丰富的图形类型和定制选项。

高度可定制：R的图形可以高度定制，从颜色、形状到布局，用户可以精细控制图形的每一个方面。

数据操作能力强：R是一种强大的数据分析工具，它可以在绘图之前对数据进行复杂的处理和分析。

集成统计分析：R语言在统计分析方面非常强大，可以直接在绘图中集成统计测试和模型结果。

动态报告：R的可视化可以轻松地集成到动态报告中，如R Markdown和Shiny，这使得结果可以交互式地呈现。

社区支持：R有一个活跃的社区，用户可以从中获得大量的教程、论坛讨论和包更新。

开源免费：R是开源软件，可以免费使用和修改，这使得它在学术界和数据科学领域非常受欢迎。

跨平台：R可以在多种操作系统上运行，包括Windows、macOS和Linux。

R语言缺点：

学习曲线：对于初学者来说，R的学习曲线可能比较陡峭，尤其是对于那些没有编程背景的用户。

性能问题：虽然R在数据处理和可视化方面表现出色，但它在处理非常大的数据集时可能会变慢。

图形更新：在R中更新图形可能需要重新运行整个脚本，这在迭代过程中可能会显得繁琐。

3D可视化有限：虽然R可以创建3D图形，但与一些专门的软件相比，它的3D可视化能力相对有限。

用户界面：R的默认用户界面（RStudio除外）可能不如一些集成开发环境（IDE）友好。

依赖管理：R的包依赖管理有时可能会导致版本冲突，需要用户手动解决。

专业图形制作：虽然R可以制作高质量的图形，但与专业的图形设计软件相比，它在制作出版级别的图形方面可能不够直观。

交互性：虽然R可以通过Shiny等工具创建交互式图形，但这些交互性可能不如专门的交互式可视化工具（如Tableau）直观和强大。

R绘图最重要的一点是无法像PS一样随意改变图形，只能在参数限定范围内修改调整。因此，要想把图做的好看，就要选择合适的绘图函数，并不停的调参。所以说，理解函数是R绘图的基础。

2.可视化

首先加载R包，创建新的文件夹，读取细胞注释后的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还可以把两个基因画在同一个图中,看右上角可以发现黄色越深的地方两个基因叠加越多：

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: 用于指定每个特征的表达值截止的向量。可以使用百分位数来指定截止点（例如，'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 : 添加显示单元格组状态的颜色栏

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包层出不穷，各种好看的图形让我们目不暇接。在不停的学习代码之余，我们还要清晰的认识到，可视化的主要目的是清晰的呈现数据的特征并揭示生物学意义。所以，我们还是要不停的看教材、看文献、学统计，不停的思考。任重而道远，我们仍需努力。下一期，我们将在此基础上，绘制饼图、堆积柱状图、箱线图、气泡图等，比较不同分组之间细胞比例差异。干货满满，欢迎大家持续追更，谢谢！