scRNA基础分析-1:安装包、导入数据、过滤质控

已经学习了转录组基本入门知识，现在了解、学习下单细胞转录组(scRNA)相关知识。
依旧是站在巨人的肩膀上看世界，入门学习目前主要参考微信公众号--生信会客厅的单细胞转录组分析教程。笔记代码，示例文件均来自推文~

scRNA基础分析-1:安装包、导入数据、过滤质控 - 简书
 scRNA基础分析-2：降维与聚类 - 简书
 scRNA基础分析-3：鉴定细胞类型 - 简书
 scRNA基础分析-4：细胞亚类再聚类、注释 - 简书
 scRNA基础分析-5：伪时间分析 - 简书
 scRNA基础分析-6：富集分析 - 简书

一、安装相关R包

install.packages('Seurat', dependencies=T)
#scRNA分析的基础包
install.packages("tidyverse", dependencies=T)
install.packages("BiocManager", dependencies=T)
BiocManager::install(c("SingleR","monocle", "DESeq2", "MAST"),ask = F,update = F) 
BiocManager::install(c("clusterProfiler","DOSE","pheatmap"),ask = F,update = F)
BiocManager::install(c("org.Hs.eg.db","org.Mm.eg.db","org.Rn.eg.db"),ask = F,update = F)

二、数据导入，创建seurat对象

library(Seurat)
scRNA.counts <- Read10X(data.dir = "filtered_feature_bc_matrix")  
#filtered_feature_bc_matrix为一个目录，里面有同一样本的三种文件。

2-1

scRNA = CreateSeuratObject(scRNA.counts)
#后续的分析大都以此为基础
table(scRNA@meta.data$orig.ident)     #查看样本的细胞数量
#1222个

三、质控，过滤不感兴趣的细胞

过滤标准有细胞基因数量、mRNA数量、核糖体基因等。因此，前期我们需要对所有细胞的相关概况进行初步探索

探索

1、细胞的基因(nFeature_RNA）、mRNA数量(nCount_RNA)

head(scRNA@meta.data)
summary(scRNA@meta.data)

3-1

2、细胞中核糖体基因比例

scRNA[["percent.mt"]] <- PercentageFeatureSet(scRNA, pattern = "^MT-")
head(scRNA@meta.data)
#计算红细胞比例，这一步没有看懂；红细胞基因比例？
HB.genes <- c("HBA1","HBA2","HBB","HBD","HBE1","HBG1","HBG2","HBM","HBQ1","HBZ")
HB_m <- match(HB.genes, rownames(scRNA@assays$RNA)) 
HB.genes <- rownames(scRNA@assays$RNA)[HB_m] 
HB.genes <- HB.genes[!is.na(HB.genes)] 
scRNA[["percent.HB"]]<-PercentageFeatureSet(scRNA, features=HB.genes)

上述可视化

#小提琴图
col.num <- length(levels(scRNA@active.ident))
# 1L
violin <- VlnPlot(scRNA,
                  features = c("nFeature_RNA", "nCount_RNA", "percent.mt"), 
                  cols =rainbow(col.num), 
                  pt.size = 0.01, #不需要显示点，可以设置pt.size = 0
                  ncol = 3)

3-2

#细胞基因数量与mRNA、核糖体基因数量的相关性
plot1 <- FeatureScatter(scRNA, feature1 = "nCount_RNA", feature2 = "nFeature_RNA")
plot2 <- FeatureScatter(scRNA, feature1 = "nCount_RNA", feature2 = "percent.mt")
plot1|plot2

3-3

质控、过滤

一般根据基因数量和核糖体比例来过滤细胞
作者建议minGene=500，maxGene=4000，pctMT=15（线粒体基因比例低于15%）。

minGene=500
maxGene=4000
pctMT=15
scRNA <- subset(scRNA, subset = nFeature_RNA > minGene & nFeature_RNA < maxGene & percent.mt < pctMT)
col.num <- length(levels(scRNA@active.ident))
violin <-VlnPlot(scRNA,
                 features = c("nFeature_RNA", "nCount_RNA", "percent.mt"), 
                 cols =rainbow(col.num), 
                 pt.size = 0.1, 
                 ncol = 3)

3-4

数据标准化与中心化（按列）

标准化是原始分数减去平均数然后除以标准差，差异基因分析会用到；

scRNA <- NormalizeData(scRNA, normalization.method = "LogNormalize", scale.factor = 10000)

中心化是原始分数减去平均数，PCA主成分分析会用到。

scale.genes <-  rownames(scRNA)
#如果自己数据集太大，也可到后面只对高变基因中心化
scRNA <- ScaleData(scRNA, features = scale.genes)

经过标准化和中心化之后，加上原始表达矩阵已经产生了三套基因表达数据

#原始表达矩阵
GetAssayData(scRNA,slot="counts",assay="RNA") 
#标准化之后的表达矩阵              
GetAssayData(scRNA,slot="data",assay="RNA")
#中心化之后的表达矩阵 
GetAssayData(scRNA,slot="scale.data",assay="RNA")

最后保存结果，供后面的分析

saveRDS(scRNA, file="scRNA.rds")

最后编辑于：2020.09.10 16:33:44