4.0 安装包
install.packages('Seurat')
加载所需要的包
library(Seurat)
library(tidyverse)
library(dplyr)
library(patchwork)
清除环境，
rm(list=ls())
读取10x的数据
BC10 文件夹内容
4640_0.png
scRNA.counts=Read10X("F:/006_data/BC10")
class(scRNA.counts)
scRNA.counts是一个S4对象。S4对象就像一颗大树一样，有主枝干，有分枝。
4642_0.png
创建Seurat对象
?CreateSeuratObject
scRNA = CreateSeuratObject(scRNA.counts ,min.cells = 3,project="os", min.features = 300)
view(scRNA)
min.cells 每个基因至少要在3个细胞中表达
min.features 每个细胞至少有多少基因表达/被检测到
查看样本的细胞数量
第一种s4对象的提取方法 点击白框
第二提取s4对象的方法 @ $交替使用
一般来说S4第一个用@
先输入枝干scRNA@assays$RNA@counts@Dim
table(scRNA@meta.data$orig.ident)
计算质控指标，去除低质量的细胞
计算细胞中线粒体基因比例，线粒体基因公认是一定要去除的。
1 线粒体是独立遗传的，不是染色体上基因控制的。
2 单细胞测序一般都是新鲜的组织，防止mRNA降解，临床新鲜的。防止线粒体快速扩增。 不新鲜，线粒体不会占比太高，线粒体占比太高，细胞可能凋亡或坏死。
scRNA[[]]
scRNA[["percent.mt"]] <- PercentageFeatureSet(scRNA, pattern = "^MT-")
x <- PercentageFeatureSet(scRNA, pattern = "^MT-")
x列名是BULK细胞标签，行名是线粒体在每个基因的占比
4646_1.png
线粒体基因都是用MT-开头命名的。用正则表达式^符号表示从首字母开始MT进行匹配提取。
计算红细胞比例，红细胞没有细胞核，没有转录组
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))
Feature、count、线粒体基因、红细胞基因占比可视化。
violin <- VlnPlot(scRNA,
features = c("nFeature_RNA", "nCount_RNA", "percent.mt","percent.HB"),
cols =rainbow(col.num),
pt.size = 0.01, #不需要显示点，可以设置pt.size = 0
ncol = 4) +
theme(axis.title.x=element_blank(), axis.text.x=element_blank(), axis.ticks.x=element_blank())
把图片画到画板上面
violin
4644_0.png
保存
ggsave("vlnplot_before_qc.pdf", plot = violin, width = 12, height = 6)
ggsave("vlnplot_before_qc.png", plot = violin, width = 12, height = 6)
这几个指标之间的相关性。 把图画到画板上，然后手动保存
plot1=FeatureScatter(scRNA, feature1 = "nCount_RNA", feature2 = "percent.mt")
plot2=FeatureScatter(scRNA, feature1 = "nCount_RNA", feature2 = "nFeature_RNA")
plot3=FeatureScatter(scRNA, feature1 = "nCount_RNA", feature2 = "percent.HB")
pearplot <- CombinePlots(plots = list(plot1, plot2, plot3), nrow=1, legend="none")
看画板
plot1
plot2
plot3
pearplot
我们可以看到，nFeature_RNA的范围在0到8000之内，每个细胞基因表达量
percent.mt代表线粒体含量
percent.HB代表红细胞含量
我们默认线粒体含量至少要小于10%，这是根据生物学知识得出的默认阈值。红细胞的数目要至少小于3%
至于nFeature_RNA和nCount_RNA的阈值怎么确定，这个要结合 pearplot的图来判断。我们质控的目标就是删除离异值。而且注意阈值尽可能取的宽松一下，防止后面分析想要的细胞得不到。
接下来从pearplot的图片来做质控---剔除离异值
nFeature_RNA选择大于300 小于7000的 nFeature_RNA选择小于100000，percent.mt小于10，percent.HB小于3
scRNA1 <- subset(scRNA, subset = nFeature_RNA > 300& nFeature_RNA < 7000 & percent.mt < 10 & percent.HB < 3 & nCount_RNA < 100000)
scRNA
scRNA1
在控制台中我们可以看到有500多细胞过滤了
过滤完之后 我们就要对数据进行均一化，使用NormalizeData这个函数。
注意均一化是用NormalizeData，标准化是用ScaleData
scRNA1 <- NormalizeData(scRNA1, normalization.method = "LogNormalize", scale.factor = 10000)
好了，这一节数据加载、质控的内容就算是做完了。
在我们关闭rstudio之前 先把环境中运行好的数据保存一下
数据将保存在之前设定好的路径中。还有保存的scRNA1，不是scRNA，因为scRNA1才是过滤好的数据。
save(scRNA1,file='scRNA1.Rdata')
