### 使用 R 语言进行单细胞数据降维聚类并实现细胞类型注释
#### 数据准备
为了执行有效的降维聚类分析，在开始之前需准备好必要的输入文件，这些通常包括 `genes.tsv` 文件（或称为 `features.tsv`）、`barcodes.tsv` 和 `matrix.mtx` 表达量矩阵[^2]。
#### 引入相关库
在 R 中操作时，Seurat 是一个广泛使用的包来处理单细胞 RNA 测序(scRNA-seq)的数据。安装 Seurat 可通过 CRAN 或 Bioconductor 安装最新版本：
```r
if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("Seurat")
```
加载所需的软件包：
```r
library(Seurat)
library(ggplot2)
```
#### 创建 Seurat 对象与初步预处理
读取表达矩阵并将之转换成适合进一步分析的形式。这一步骤涉及创建一个新的 Seurat 对象，并对其进行基本的质量控制(QC)，比如去除低质量的细胞和特征。
```r
# 假设已经存在名为 pbmc_small 的对象
pbmc <- CreateSeuratObject(counts = pbmc_small, project = "PBMC_3k", min.cells = 3, min.features = 200)
# 执行 QC 步骤
pbmc[["percent.mt"]] <- PercentageFeatureSet(pbmc, pattern = "^MT-")
VlnPlot(pbmc, features = c("nFeature_RNA", "nCount_RNA", "percent.mt"), ncol = 3)
FilterPlot(pbmc, idents = NULL, slot = "data") + NoLegend()
pbmc <- subset(pbmc, subset = nFeature_RNA > 200 & nCount_RNA < 5000 & percent.mt < 10)
```
#### 数据标准化、缩放以及识别高度可变基因(HVGs)
接下来是对经过筛选后的数据集应用标准化方法，之后计算每个样本中每条基因的标准差得分，从而挑选出那些表现出显著变化水平的 HVGs 来用于后续分析。
```r
pbmc <- NormalizeData(pbmc, normalization.method = "LogNormalize", scale.factor = 10000)
pbmc <- FindVariableFeatures(pbmc, selection.method = "vst", nfeatures = 2000)
top10 <- head(x = VariableFeatures(object = pbmc), 10L)
```
#### PCA 主成分分析 (PCA)
利用主成分分析技术减少维度数量的同时保留尽可能多的信息。此阶段会生成一组新的坐标轴——即所谓的“主成分”，它们能够捕捉原始变量间的最大方差。
```r
all.genes <- rownames(x = pbmc@assays$RNA@counts)
pbmc <- ScaleData(pbmc, vars.to.regress = "percent.mt")
pbmc <- RunPCA(pbmc, npcs = 30, verbose = FALSE)
Print(x = pbmc[["pca"]], dims = 1:5, reduction.key = "")
DimHeatmap(pbmc, dims = 1, cells.use = 500, balanced = TRUE, group.by = "stim")
ElbowPlot(pbmc)
```
#### 构建邻接图与 Louvain 社区检测算法构建集群
基于 k-nearest neighbors(KNN) 图模型建立细胞之间的连接关系；随后采用改进版的 Louvain 方法来进行社区发现，以此定义不同的细胞群体。
```r
pbmc <- FindNeighbors(pbmc, dims = 1:10)
pbmc <- FindClusters(pbmc, resolution = 0.5)
table(Idents(object = pbmc))
DimPlot(pbmc, reduction = "umap", label = TRUE, pt.size = 0.5)
```
#### 差异表达分析及标记基因查找
最后一步是找出区分各个簇之间差异性的标志物，这对于理解各组内的生物学特性至关重要。
```r
markers <- FindAllMarkers(pbmc, only.pos = TRUE, min.pct = 0.25, logfc.threshold = 0.25)
head(markers, n = 5)
DefaultAssay(pbmc) <- "RNA"
VlnPlot(pbmc, features = c("MS4A1", "CD79A", "MALAT1"))
DoHeatmap(pbmc, features = c("MS4A1", "CD79A", "MALAT1"))
```
对于未被成功鉴定出来的某些细胞群，可能是由于来自不同批次或者物种间存在的固有差异所造成的[^3]。
最低0.47元/天 解锁文章
扫一扫
1442
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
