单细胞Seurat流程与步骤详解

开始之前我们先看一下标准流程

# 加载必要的库
library(Seurat)
library(harmony)

# 数据预处理与分析
subset_data <- subset_data %>%
  Seurat::NormalizeData(verbose = FALSE) %>%  # 归一化数据
  FindVariableFeatures(selection.method = "vst", nfeatures = 2000) %>%  # 选择2000个变异基因
  ScaleData(verbose = FALSE) %>%  # 数据缩放
  RunPCA(npcs = 50, verbose = FALSE)  # 主成分分析（PCA）

# Harmony批次效应校正
subset_data <- subset_data %>% RunHarmony("stim", plot_convergence = TRUE)

# 获取Harmony嵌入数据
harmony_embeddings <- Embeddings(subset_data, 'harmony')

# 降维和可视化
dims = 1:30
subset_data <- subset_data %>%
  RunUMAP(reduction = "harmony", dims = dims) %>%  # 使用Harmony嵌入做UMAP
  RunTSNE(reduction = "harmony", dims = dims)  %>% # 使用Harmony嵌入做t-SNE
  FindNeighbors(reduction = "harmony", dims = dims) %>%
  Findclusters(reduction = "harmony", dims = dims)
# 可视化结果
DimPlot(subset_data, reduction = "umap")  # UMAP可视化
DimPlot(subset_data, reduction = "tsne")  # t-SNE可视化

# 如果你已经定义了细胞类型，可以用以下方式展示不同类型的细胞
DimPlot(subset_data, group.by = "cell_type")  # 基于已有的细胞类型标签绘制图

详细解释

1. 数据归一化 (NormalizeData)

subset_data = subset_data %>%
  Seurat::NormalizeData(verbose = FALSE)

• 目的：数据归一化的目的是将每个细胞的基因表达值调整到同一尺度，消除测序深度（library size）对表达数据的影响。通常，基因表达数据会进行对数变换，以便更好地比较各个基因的表达水平。
• 使用的输入数据：它使用的是 subset_data 对象中的原始基因表达矩阵（通常是原始的 RNA counts）。
• 结果：NormalizeData 函数会创建一个新的归一化数据集，结果存储在 subset_data 中。在 verbose = FALSE 时，函数不会输出过多的运行信息。

2. 选择可变基因 (FindVariableFeatures)

subset_data = subset_data %>%
  FindVariableFeatures(selection.method = "vst", nfeatures = 2000)

• 目的：选择在不同细胞之间变异性最大的基因，通常选择前 2000 个变异性最大的基因。这些基因常被用于后续分析，如降维和聚类。
• 使用的输入数据：它使用归一化后的数据（即上一步 NormalizeData 处理过的数据），通过计算每个基因在所有细胞中的变异性来选择最具代表性的基因。
• 结果：FindVariableFeatures 通过计算基因的变异性，选出变异性最大的 2000 个基因，并将这些基因的信息存储在 subset_data 的 @meta.data 中。

3. 数据缩放 (ScaleData)

subset_data = subset_data %>%
  ScaleData(verbose = FALSE)

• 目的：将选择的变异基因数据进行标准化，使得每个基因的均值为 0，标准差为 1。这样可以确保每个基因在降维（如 PCA）时具有相同的贡献，不会因为基因的表达水平差异过大而影响结果。
• 使用的输入数据：它使用的是上一步选择的变异基因数据，通常是经过归一化的数据。
• 结果：ScaleData 处理后，每个基因的数据被缩放，结果存储在 subset_data 对象中。

4. 主成分分析 (PCA) (RunPCA)

subset_data = subset_data %>%
  RunPCA(npcs = 50, verbose = FALSE)

• 目的：主成分分析（PCA）是降维的一种方法，它通过计算数据中的主成分（PCs）来减少特征空间的维度。在这里，我们选择了前 50 个主成分（npcs = 50）。
• 使用的输入数据：它使用的是缩放后的数据（即 ScaleData 处理过的数据）。PCA 会基于这些数据计算出每个细胞在每个主成分上的投影。
• 结果：PCA 结果被存储在 subset_data 对象中，作为一个新的维度数据，可以用于后续的降维和可视化步骤。

5. 绘制肘部图 (ElbowPlot)

ElbowPlot(subset_data, ndims = 50)

• 目的：肘部图用来帮助选择适当的主成分数目，通常选择肘部（即误差变化急剧减缓的地方）之前的主成分数量作为分析的维度数。
• 使用的输入数据：它使用的是前一步的 PCA 结果（即 RunPCA 生成的主成分数据）。
• 结果：肘部图会展示不同主成分数与方差解释度的关系，帮助你判断哪些主成分对数据变异贡献较大，通常选择那些方差解释度较高的主成分。

6. 绘制小提琴图 (VlnPlot)

VlnPlot(subset_data, features = c("nFeature_RNA", "nCount_RNA", "percent.mt"), ncol = 3)

• 目的：通过小提琴图（violin plot）展示不同细胞的基因数目（nFeature_RNA）、转录本数目（nCount_RNA）以及线粒体基因的百分比（percent.mt）。这些指标常用于质量控制，帮助识别低质量细胞或过度激活的细胞。
• 使用的输入数据：它使用 subset_data 对象中的元数据（@meta.data）以及基因表达数据。
• 结果：图形会显示这些不同特征在细胞中的分布情况，帮助你了解数据的质量。

7. 数据整合（Harmony）

library('harmony')
subset_data <- subset_data %>%
  RunHarmony("stim", plot_convergence = TRUE)

• 目的：使用 Harmony 方法整合批次效应（batch effects）或条件效应（如刺激 vs. 对照），使得整合后的数据不受这些效应的干扰。stim 是数据中的一个元数据列，用来标记不同的实验条件（如实验组 vs. 对照组）。
• 使用的输入数据：它默认使用的是 subset_data 对象中的 PCA 结果（即 RunPCA 或其他降维方法的结果）。
• 结果：RunHarmony 会基于所指定的元数据（在这里是 stim）来消除批次效应或实验条件效应，生成新的“和谐”嵌入（harmony embeddings），它们存储在 subset_data 对象中。

8. 降维与可视化（UMAP 和 t-SNE）

subset_data <- subset_data %>%
  RunUMAP(reduction = "harmony", dims = dims) %>%
  RunTSNE(reduction = "harmony", dims = dims) %>%
  FindNeighbors(reduction = "harmony", dims = dims)

• 目的：使用 UMAP 和 t-SNE 进行非线性降维可视化。UMAP 和 t-SNE 都是常用于高维数据可视化的技术，通过将数据从高维空间降到 2D 或 3D 空间，便于观察细胞群体的分布。FindNeighbors 则用于根据指定的维度（如前 30 个 Harmony 维度）计算细胞之间的相似性，生成邻居图。
• 使用的输入数据：RunUMAP、 RunTSNE使用的是 Harmony 处理过的数据（即 RunHarmony 的结果），并基于指定的维度（dims）进行降维。
• 结果：生成的 UMAP 和 t-SNE 图将存储在 subset_data 对象中，帮助你直观地查看细胞在低维空间的分布。

9. 细胞聚类 (FindClusters)

subset_data = FindClusters(subset_data, resolution = 0.3)

• 目的：进行细胞聚类。聚类方法（如 Louvain 或 Leiden）将细胞分组，目的是识别具有相似基因表达模式的细胞群体。
• 使用的输入数据：它使用的是 FindNeighbors 生成的细胞邻接图以及降维后的数据。
• 结果：细胞被分配到不同的群体（clusters），并存储在 subset_data@meta.data 中的 seurat_clusters 列。

注意：

FindNeighbors 和 FindClusters 主要是用来进行细胞聚类的，它们会基于数据的相似性计算细胞之间的邻近关系，并通过聚类算法（如 Louvain 或 Leiden 方法）将细胞分成不同的群体。这对于数据的探索和分析非常重要，尤其是当你没有预定义细胞类型时。

但是，如果你已经定义好了细胞类型，并且不打算根据数据的相似性进行自动聚类，那么确实不需要执行 FindNeighbors 和 FindClusters。在这种情况下，你可以跳过这两步，只使用你定义好的细胞类型标签来进行后续的分析和可视化。

• FindNeighbors：计算细胞之间的相似性（通常使用主成分分析（PCA）结果）来为后续的聚类步骤提供邻近关系。如果你已经有了细胞类型标签，这一步不是必需的。
• FindClusters：基于 FindNeighbors 计算的邻近关系进行细胞聚类，将相似的细胞分为同一类。如果你已经定义了细胞类型，则不需要根据数据进行重新聚类。

但如果你想使用你已经定义好的细胞类型进行可视化（如 UMAP），可以直接使用这些标签进行颜色映射，无需运行聚类过程。只需在 UMAP 或其他可视化方法中使用细胞类型标签来显示不同的群体。

举个例子：

假设你已经有了细胞类型信息，并且想直接画出基于这些细胞类型的 UMAP 图：

# 假设细胞类型已存储在 subset_data@meta.data$cell_type 中
UMAP_plot <- DimPlot(subset_data, group.by = "cell_type")

这样，你就可以在 UMAP 图中看到基于已知细胞类型的分布，而无需运行 FindNeighbors 和 FindClusters。

总之，如果你已经有了细胞类型的信息，并且不需要基于数据进行新的聚类，FindNeighbors 和 FindClusters 这两步可以省略。

总结

本文代码块完成了从数据归一化、选择变异基因、PCA降维到批次效应校正和细胞聚类、ump和tsne的完整流程。每一步都在前一步的结果基础上进行，确保数据经过预处理和整合，最终得到可用于细胞分群、降维可视化的高质量数据。