刘小泽写于2020.7.21
为何取名叫“交响乐”？因为单细胞分析就像一个大乐团，需要各个流程的协同配合
单细胞交响乐1-常用的数据结构SingleCellExperiment
单细胞交响乐2-scRNAseq从实验到下游简介
单细胞交响乐3-细胞质控
单细胞交响乐4-归一化
单细胞交响乐5-挑选高变化基因
单细胞交响乐6-降维
单细胞交响乐7-聚类分群
单细胞交响乐8-marker基因检测
单细胞交响乐9-细胞类型注释
单细胞交响乐9-细胞类型注释
单细胞交响乐10-数据集整合后的批次矫正
单细胞交响乐11-多样本间差异分析
单细胞交响乐12-检测Doublet
单细胞交响乐13-细胞周期推断
单细胞交响乐14-细胞轨迹推断
单细胞交响乐15-scRNA与蛋白丰度信息结合
单细胞交响乐16-处理大型数据
单细胞交响乐17-不同单细胞R包的数据格式相互转换
单细胞交响乐18-实战一 Smart-seq2
单细胞交响乐19-实战二 STRT-Seq
单细胞交响乐20-实战三 10X 未过滤的PBMC数据
单细胞交响乐21-实战三 批量处理并整合多个10X PBMC数据
单细胞交响乐22-实战五 CEL-seq2
单细胞交响乐23-实战六 CEL-seq
单细胞交响乐24-实战七 SMARTer 胰腺细胞
单细胞交响乐25-实战八 Smart-seq2 胰腺细胞
单细胞交响乐26-实战九 胰腺细胞数据整合
单细胞交响乐27-实战十 CEL-seq-小鼠造血干细胞
单细胞交响乐28-实战十一 Smart-seq2-小鼠造血干细胞
1 前言
前面的种种都是作为知识储备，但是不实战还是记不住前面的知识
这是第十二个实战练习
数据来自Pijuan-Sala et al. (2019)，研究的是小鼠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
)
欢迎关注我们的公众号～_～
我们是两个农转生信的小硕，打造生信星球，想让它成为一个不拽术语、通俗易懂的生信知识平台。需要帮助或提出意见请后台留言或发送邮件到jieandze1314@gmail.com
