一直在路上,人生才会充满无限可能。
先定个小目标,系统的学习完刘小泽老师的《单细胞交响乐1-31》
两张非常重要的图
原来不用Seurat
也可以做到质控、挑高变化基因:
下面分析下droplet-based的视网膜数据【Macosko et al. (2015)】
示例数据
rm(list = ls())
# options("repos"= c(CRAN="https://mirrors.tuna.tsinghua.edu.cn/CRAN/"))
# options(BioC_mirror="http://mirrors.tuna.tsinghua.edu.cn/bioconductor/")
library(scRNAseq)
library(scater)
library(scran)
# sce <- MacoskoRetinaData()
# save(sce,file = "MacoskoRetinaData.RData")
load("MacoskoRetinaData.RData")
质控和归一化
# 质控
is.mito <- grepl("^MT-",rownames(sce)) #去除线粒体基因
## 去除低质量基因
qcstats <- perCellQCMetrics(sce,subsets=list(Mito=is.mito))
filtered <- quickPerCellQC(qcstats,percent_subsets="subsets_Mito_percent") #去除
sce <- sce[,!filtered$discard]
# 归一化
sce <- logNormCounts(sce)
sce
找高变基因
# 找高变基因
dec <- scran::modelGeneVar(sce)
hvg <- scran::getTopHVGs(dec,prop = 0.1)
dec
hvg[1:10];dim(hvg)
降维聚类
# 降维
set.seed(1234)
sce <- runPCA(sce,ncomponents=25,subset_row=hvg)
sce <- runUMAP(sce,dimred="PCA",external_neighbors=T)
#聚类
g <- buildSNNGraph(sce,use.dimred="PCA")
sce$cluster <- factor(igraph::cluster_louvain(g)$membership)
可视化
plotUMAP(sce, colour_by="cluster")
下面就该一套下游分析了,差异分析、亚群注释、细胞周期推断、细胞互做等。
参考链接: