1 前言
数据来自，研究的是小鼠E8.5发育阶段的嵌合胚胎
数据准备
# 自己下载
library(MouseGastrulationData)
sce.chimera <- WTChimeraData(samples=5:10)
# 或者加载之前分享的数据
load('sce.chimera.RData')
sce.chimera
## class: SingleCellExperiment
## dim: 29453 20935
## metadata(0):
## assays(1): counts
## rownames(29453): ENSMUSG00000051951 ENSMUSG00000089699 ...
## ENSMUSG00000095742 tomato-td
## rowData names(2): ENSEMBL SYMBOL
## colnames(20935): cell_9769 cell_9770 ... cell_30702 cell_30703
## colData names(11): cell barcode ... doub.density sizeFactor
## reducedDimNames(2): pca.corrected.E7.5 pca.corrected.E8.5
## altExpNames(0):
names(colData(sce.chimera))
# [1] "cell" "barcode"
# [3] "sample" "stage"
# [5] "tomato" "pool"
# [7] "stage.mapped" "celltype.mapped"
# [9] "closest.cell" "doub.density"
简单看一下colData中的各个信息
其中包含了6个样本的信息，总共20935个细胞
table(sce.chimera$sample)
#
# 5 6 7 8 9 10
# 2298 1026 2740 2904 4057 6401
整合行名
library(scater)
rownames(sce.chimera) <- uniquifyFeatureNames(
rowData(sce.chimera)$ENSEMBL, rowData(sce.chimera)$SYMBOL)
2 简单质控
之前作者已经对数据进行了质控，并把细胞做了标志，这里只需要把标记“stripped”、“Doublet”的细胞去掉即可
drop <- sce.chimera$celltype.mapped %in% c("stripped", "Doublet")
table(drop)
# drop
# FALSE TRUE
# 19426 1509
sce.chimera <- sce.chimera[,!drop]
3 归一化
看到原来数据中也计算了size factors，那么这里就不需要计算，直接应用
sce.chimera <- logNormCounts(sce.chimera)
4 找表达量高变化基因
我们的数据有6个样本，可以说异质性非常高了。把它们当做不同的批次信息，并尽可能多地从中保存基因
library(scran)
dec.chimera <- modelGeneVar(sce.chimera, block=sce.chimera$sample)
chosen.hvgs <- dec.chimera$bio > 0
table(chosen.hvgs)
# chosen.hvgs
# FALSE TRUE
# 14754 14699
5 数据整合并矫正批次效应
使用了一种“层次整合”的方法，就是先将同种表型样本整合起来（比如3个处理和3个对照先内部整合），再将不同表型的样本组合（将处理和对照整合）
correctExperiments的含义是：Apply a correction to multiple SingleCellExperiment objects,
library(batchelor)
set.seed(01001001)
# 下面的merge.order就设置了整合的顺序
merged <- correctExperiments(sce.chimera,
batch=sce.chimera$sample,
subset.row=chosen.hvgs,
PARAM=FastMnnParam(
merge.order=list(
list(1,3,5), # WT (3 replicates)
list(2,4,6) # td-Tomato (3 replicates)
)
)
)
看下结果：lost.var 值越大表示丢失的真实生物异质性越多
metadata(merged)$merge.info$lost.var
## 5 6 7 8 9 10
## [1,] 0.000e+00 0.0204433 0.000e+00 0.0169567 0.000000 0.000000
## [2,] 0.000e+00 0.0007389 0.000e+00 0.0004409 0.000000 0.015474
## [3,] 3.090e-02 0.0000000 2.012e-02 0.0000000 0.000000 0.000000
## [4,] 9.024e-05 0.0000000 8.272e-05 0.0000000 0.018047 0.000000
## [5,] 4.321e-03 0.0072518 4.124e-03 0.0078280 0.003831 0.007786
Large proportions of lost variance (>10%) suggest that correction is removing genuine biological heterogeneity.
6 聚类
g <- buildSNNGraph(merged, use.dimred="corrected")
clusters <- igraph::cluster_louvain(g)
colLabels(merged) <- factor(clusters$membership)
看分群与细胞类型之间关系
tab <- table(Cluster=colLabels(merged), Sample=merged$sample)
library(pheatmap)
pheatmap(log10(tab+10), color=viridis::viridis(100))
7 降维
merged <- runTSNE(merged, dimred="corrected")
merged <- runUMAP(merged, dimred="corrected")
gridExtra::grid.arrange(
plotTSNE(merged, colour_by="label", text_by="label", text_col="red"),
plotTSNE(merged, colour_by="batch"),
ncol=2
)
