跟着小鱼头学单细胞测序-单细胞数据的整合

# Seurat-merge示例

# 示例1：合并两个10X PBMC数据集
# 首先读入并创建两个seurat objects: pbmc4k 和pbmc8k

library(Seurat)
pbmc4k.data <- Read10X(data.dir = "../pbmc4k/filtered_gene_bc_matrices/GRCh38/")
pbmc4k <- CreateSeuratObject(counts = pbmc4k.data, project = "PBMC4K")
pbmc4k

pbmc8k.data <- Read10X(data.dir = "../pbmc8k/filtered_gene_bc_matrices/GRCh38/")
pbmc8k <- CreateSeuratObject(counts = pbmc8k.data, project = "PBMC8K")
pbmc8k

# 合并
# merge()函数会合并两个对象的原始表达数据；两个数据集的细胞barcode可能会重复，
# 为了区分细胞来源，通过`add.cell.ids`给不同样本的细胞barcode添加一个id

pbmc.combined <- merge(pbmc4k, y = pbmc8k, add.cell.ids = c("4K", "8K"), project =           "PBMC12K")
pbmc.combined
head(colnames(pbmc.combined))
table(pbmc.combined$orig.ident)

# 示例2：合并多个数据集：pbmc3k pbmc4k pbmc8k

InstallData("pbmc3k")
pbmc3k <- LoadData("pbmc3k", type = "pbmc3k.final")
pbmc3k
pbmc.big <- merge(pbmc3k, y = c(pbmc4k, pbmc8k), add.cell.ids = c("3K", "4K", "8K"), project = "PBMC15K")
pbmc.big
table(pbmc.big$orig.ident)

#end.

# CCA+MNN 整合

# 示例:
library(Seurat)
library(SeuratData)
library(patchwork)

# 加载数据
InstallData("ifnb")
LoadData("ifnb")

# 将数据分为两个seurat 对象 ：stim and CTRL
ifnb.list <- SplitObject(ifnb, split.by = "stim")

# 分别对两个数据集进行标准化并找出高变基因
ifnb.list <- lapply(X = ifnb.list, FUN = function(x) {
    x <- NormalizeData(x)
    x <- FindVariableFeatures(x, selection.method = "vst", nfeatures = 2000)
})

# 选择在数据集中重复的高变基因用于整合
features <- SelectIntegrationFeatures(object.list = ifnb.list)

#整合：基于高变基因找出两两数据集间的“锚点细胞”
immune.anchors <- FindIntegrationAnchors(object.list = ifnb.list, anchor.features = features)
# 生成一个新的 'integrated' assay
immune.combined <- IntegrateData(anchorset = immune.anchors)

# 通过 'integrated' assay的整合数据进行下游的细胞类型分类等分析；原始的表达数据仍保存在 'RNA' assay
DefaultAssay(immune.combined) <- "integrated"
immune.combined <- ScaleData(immune.combined, verbose = FALSE)
immune.combined <- RunPCA(immune.combined, npcs = 30, verbose = FALSE)
immune.combined <- RunUMAP(immune.combined, reduction = "pca", dims = 1:30)
immune.combined <- FindNeighbors(immune.combined, reduction = "pca", dims = 1:30)
immune.combined <- FindClusters(immune.combined, resolution = 0.5)

# 可视化：查看两个数据集之间细胞分类情况
DimPlot(immune.combined, reduction = "umap", split.by = "stim")

#end.

# 示例

library(Seurat)

# 分别对每个数据集进行进行标准化并找出高变基因
bm280k.data <- Read10X_h5("../data/ica_bone_marrow_h5.h5")
bm280k <- CreateSeuratObject(counts = bm280k.data, min.cells = 100, min.features = 500)
bm280k.list <- SplitObject(bm280k, split.by = "orig.ident")
bm280k.list <- lapply(X = bm280k.list, FUN = function(x) {
    x <- NormalizeData(x, verbose = FALSE)
    x <- FindVariableFeatures(x, verbose = FALSE)
})

# 选择在数据集中重复的高变基因用于整合
features <- SelectIntegrationFeatures(object.list = bm280k.list)

# 分别对每个数据集进行PCA降维，用于下面的rpca方法
bm280k.list <- lapply(X = bm280k.list, FUN = function(x) {
    x <- ScaleData(x, features = features, verbose = FALSE)
    x <- RunPCA(x, features = features, verbose = FALSE)
})

# 利用reference 和rpca寻找锚点细胞，进行整合
anchors <- FindIntegrationAnchors(object.list = bm280k.list, reference = c(1, 2), reduction = "rpca", dims = 1:50)
bm280k.integrated <- IntegrateData(anchorset = anchors, dims = 1:50)

# 整合后，查看数据集的分布状况
bm280k.integrated <- ScaleData(bm280k.integrated, verbose = FALSE)
bm280k.integrated <- RunPCA(bm280k.integrated, verbose = FALSE)
bm280k.integrated <- RunUMAP(bm280k.integrated, dims = 1:50)
DimPlot(bm280k.integrated, group.by = "orig.ident")

# end.