scRNA-seq单细胞学习单细胞分析专题

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

2020-08-25  本文已影响0人  小贝学生信

已经学习了转录组基本入门知识,现在了解、学习下单细胞转录组(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个

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

探索

head(scRNA@meta.data)
summary(scRNA@meta.data)
3-1
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
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)
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")
上一篇下一篇

猜你喜欢

热点阅读