2024-06-16

探索细胞的奇妙之旅：splatter带你轻松玩转单细胞测序数据！（二）

原创小果生信果 2023-10-12 19:00:17

收录于话题

在上一篇教程中，小果教大家如何使用splatter和scater包进行单细胞测序数据的模拟和分析，其中包括模拟数据的快速生成、参数估计、自定义参数创建、数据提取、数据标准化和降维处理等，那么今天小果来继续带大家学习这个R包的使用方法，今天的学习也将基于昨天的代码哦，想要深入学习的小伙伴和小果一起看下去吧！

首先，我们来回顾一下上一篇中使用splatter包的代码：

library(splatter)library(scater)set.seed(1)sce <- mockSCE() # class: SingleCellExperimentcounts(sce)params <- splatEstimate(sce)  # Params objectparams <- newSplatParams() #默认10000 Genes，100 CellsgetParam(params, "nGenes")getParam(params, "nCells")getParams(params, c("nGenes", "mean.rate", "mean.shape"))params <- setParam(params, "nGenes", 5000)params <- setParam(params, "batchCells", 50)params <- setParams(params, mean.shape = 0.5, de.prob = 0.2)params <- setParams(params, update = list(nGenes = 8000, mean.rate = 0.5))sim <- splatSimulate(params)class(sim)head(rowData(sim))head(colData(sim))names(assays(sim))assays(sim)$CellMeans[1:5, 1:5]counts <- counts(sim)counts[1:3,1:5]class(counts) # "matrix" "array" typeof(counts)  # [1] "integer"dim(counts) #[1] 8000  100sim <- logNormCounts(sim)counts(sim)[1:3,1:5]logcounts(sim)[1:3,1:5]sim <- runPCA(sim)plotPCA(sim)sim.groups <- splatSimulate(group.prob = c(0.3, 0.7), method = "groups",                            verbose = FALSE)sim.groups <- logNormCounts(sim.groups)sim.groups <- runPCA(sim.groups)plotPCA(sim.groups, colour_by = "Group")
dim(counts(sim.groups))rowData(sim.groups)$DEFacGroup1rowData(sim.groups)$DEFacGroup2metadata(sim.groups)sim.paths <- splatSimulate(de.prob = 0.2, nGenes = 1000, method = "paths",                           verbose = FALSE)sim.paths <- logNormCounts(sim.paths)sim.paths <- runPCA(sim.paths)plotPCA(sim.paths, colour_by = "Step")colData(sim.paths)$Stepsim.batches <- splatSimulate(batchCells = c(50, 50), verbose = FALSE)sim.batches <- logNormCounts(sim.batches)sim.batches <- runPCA(sim.batches)plotPCA(sim.batches, colour_by = "Batch")rowData(sim.batches)$BatchFacBatch1rowData(sim.batches)$BatchFacBatch2dev.off()

继续上面的操作，我们来继续学习这个包的应用：

比较SingleCellExperiment对象

首先我们使用splatSimulate函数生成两个SingleCellExperiment对象sim1和sim2，分别表示两个不同的实验条件。然后使用compareSCEs函数将这两个对象进行比较，并将结果保存在comparison变量中。

### 7.比较SingleCellExperiment对象## 相互比较sim1 <- splatSimulate(nGenes = 1000, batchCells = 20, verbose = FALSE)sim2 <- simpleSimulate(nGenes = 1000, nCells = 20, verbose = FALSE)comparison <- compareSCEs(list(Splat = sim1, Simple = sim2))

接着，我们来学习如何访问comparison对象中的不同属性，包括RowData和ColData，以及Plots中的Means、LibrarySizes和ZerosCell。

names(comparison)comparison$RowDatacomparison$ColData
names(comparison$Plots)comparison$Plots$Meanscomparison$Plots$LibrarySizescomparison$Plots$ZerosCell

之后，代码再次生成三个SingleCellExperiment对象，分别命名为sim1、sim2和sim3。然后使用diffSCEs函数将这三个对象与参考对象“Simple”进行比较，并将结果保存在difference变量中。

## 与参考进行比较sim1 <- splatSimulate(nGenes = 1000, batchCells = 100, verbose = FALSE)sim2 <- splatSimulate(nGenes = 1000, batchCells = c(40, 60), verbose = FALSE)sim3 <- simpleSimulate(nGenes = 1000, nCells = 100, verbose = FALSE)difference <- diffSCEs(list(Splat1 = sim1, Splat2 = sim2, Simple = sim3), ref = "Simple")difference$Plots$Meansdifference$QQPlots$Means

批次效应参数设置

接下来，我们学习如何使用splatSimulate函数生成带有批次效应的数据。以下代码分别生成了两个SingleCellExperiment对象sim1和sim2，使用不同的批次效应参数。然后对这两个对象进行了log2标准化和PCA分析，并使用plotPCA函数将结果可视化。

library("splatter")library("scater")library("ggplot2")
# Simulation with small batch effectssim1 <- splatSimulate(params, batchCells = c(100, 100),                      batch.facLoc = 0.001, batch.facScale = 0.001,                      verbose = FALSE)sim1 <- logNormCounts(sim1)sim1 <- runPCA(sim1)plotPCA(sim1, colour_by = "Batch") + ggtitle("Small batch effects")
# Simulation with big batch effectssim2 <- splatSimulate(params, batchCells = c(100, 100),                      batch.facLoc = 0.5, batch.facScale = 0.5,                      verbose = FALSE)sim2 <- logNormCounts(sim2)sim2 <- runPCA(sim2)plotPCA(sim2, colour_by = "Batch") + ggtitle("Big batch effects")

我们分别来看下两次可视化后的结果吧！

接着，小果向大家展示了如何设置批次效应参数来消除批次效应。首先生成一个带有批次效应的数据sim1，然后不移除批次效应进行log2标准化和PCA分析；接着生成一个没有批次效应的数据sim2，然后对其进行log2标准化和PCA分析。最后，使用plotPCA函数将结果可视化并进行对比。

# 是否消除批次效应设置sim1 <- splatSimulate(params, batchCells = c(100, 100), batch.rmEffect = FALSE,                      verbose = FALSE)sim1 <- logNormCounts(sim1)sim1 <- runPCA(sim1)plotPCA(sim1, colour_by = "Batch") + ggtitle("With batch effects")
sim2 <- splatSimulate(params, batchCells = c(100, 100), batch.rmEffect = TRUE,                      verbose = FALSE)sim2 <- logNormCounts(sim2)sim2 <- runPCA(sim2)plotPCA(sim2, colour_by = "Batch") + ggtitle("Batch effects removed")

可视化后的结果如下：

离群值参数设置

接下来，小果将展示了如何设置离群值参数，并使用splatSimulate函数生成带有离群值的数据。以下代码分别生成了两个SingleCellExperiment对象sim1和sim2，使用不同的离群值概率参数。

### 9. 离群值参数设置## ----outlier-prob-------------------------------------------------------------# Few outlierssim1 <- splatSimulate(out.prob = 0.001, verbose = FALSE)ggplot(as.data.frame(rowData(sim1)),       aes(x = log10(GeneMean), fill = OutlierFactor != 1)) +  geom_histogram(bins = 100) +  ggtitle("Few outliers")
# Lots of outlierssim2 <- splatSimulate(out.prob = 0.2, verbose = FALSE)ggplot(as.data.frame(rowData(sim2)),       aes(x = log10(GeneMean), fill = OutlierFactor != 1)) +  geom_histogram(bins = 100) +  ggtitle("Lots of outliers")

最后绘图结果让我们来一起看一下：

设置组参数

接下来，小果将教大家如何设置组参数，并使用splatSimulateGroups函数生成具有不同组的数据。首先我们恶意生成一个包含一个小组和一个大组的数据sim1，然后对其进行log2标准化和PCA分析；接着生成一个包含五个组的数据sim2，然后对其进行log2标准化和PCA分析。最后，使用plotPCA函数将结果可视化并进行对比。

### 10.设置组参数# One small group, one big groupparams.groups <- newSplatParams(batchCells = 500, nGenes = 1000)sim1 <- splatSimulateGroups(params.groups, group.prob = c(0.9, 0.1),                            verbose = FALSE)sim1 <- logNormCounts(sim1)sim1 <- runPCA(sim1)plotPCA(sim1, colour_by = "Group") + ggtitle("One small group, one big group")
# Five groupssim2 <- splatSimulateGroups(params.groups,                            group.prob = c(0.2, 0.2, 0.2, 0.2, 0.2),                            verbose = FALSE)sim2 <- logNormCounts(sim2)sim2 <- runPCA(sim2)plotPCA(sim2, colour_by = "Group") + ggtitle("Five groups")

可视化后的结果如下：

最后，小果将展示了如何设置差异表达基因参数，并使用splatSimulateGroups函数生成具有不同差异表达基因的数据。

### 11. 设置差异表达基因参数## de.prob 参数# Few DE genesparams.groups <- newSplatParams(batchCells = 500, nGenes = 1000)sim1 <- splatSimulateGroups(params.groups, group.prob = c(0.5, 0.5),                            de.prob = 0.01, verbose = FALSE)sim1 <- logNormCounts(sim1)sim1 <- runPCA(sim1)plotPCA(sim1, colour_by = "Group") + ggtitle("Few DE genes")
# Lots of DE genessim2 <- splatSimulateGroups(params.groups, group.prob = c(0.5, 0.5),                            de.prob = 0.3, verbose = FALSE)sim2 <- logNormCounts(sim2)sim2 <- runPCA(sim2)plotPCA(sim2, colour_by = "Group") + ggtitle("Lots of DE genes")