上一篇说道SCDE包进行单细胞测序的差异表达分析。
这一篇就来说说DEsingle这个包。
1、安装
先上包的官网链接
https://github.com/miaozhun/DEsingle/
安装方法
if(!require(BiocManager)) install.packages("BiocManager")
BiocManager::install("DEsingle")
2、示例数据预分析
载入包和数据
library(DEsingle)
data(TestData)
dim(counts)
counts[1:6, 1:6]
length(group)
summary(group)
分组。
划重点：这个包奔溃的一点是只能比较两个组。如果srurat出来还几个组，那就只能两两比较。
或者，如果各位大神有其他妙招，欢迎留言交流学习。
# Specifying the two groups to be compared
# The sample number in group 1 and group 2 is 50 and 100 respectively
group <- factor(c(rep(1,50), rep(2,100)))
# Detecting the DE genes
results <- DEsingle(counts = counts, group = group)
# Dividing the DE genes into 3 categories at threshold of FDR < 0.05
results.classified <- DEtype(results = results, threshold = 0.05)
3、个人的具体数据分析
自己的数据是10×genomics产生的，所以势必要到Seurat包里面分析一遍。
但是经过Seurat包处理的数据是Seurat格式。
而DEsingle这个包只可以识别SingleCellExperiment格式的数据。
那就。。。。哎~~转换格式吧。
转换格式很简单。
借助于library(SingleCellExperiment)，一键转换。
library(SingleCellExperiment)自己去搜索下吧，不难。
接下来，导入自己的数据。
library(Seurat)
library(SingleCellExperiment)
library(DEsingle)
pbmc.data <- Read10X(data.dir = "F:/Seurat/")
pbmc <- CreateSeuratObject(counts = pbmc.data, project = "pbmc", min.cells = 3, min.features = 200)
pbmc
### 预处理与质量控制
#计算线粒体RNA比例
pbmc[["percent.mt"]] <- PercentageFeatureSet(pbmc, pattern = "^MT-")
#去除表达基因数量少于500或多于6000的细胞，去除线粒体RNA比例高于40%的细胞.这个具体标准的筛选，根据之前violin plot来确定。
pbmc <- subset(pbmc, subset = nFeature_RNA > 2000 & nFeature_RNA < 6000 & percent.mt < 30)
dim(pbmc)
##或者直接读入之前保存的数据
pbmc <- readRDS(file = "F:/Seurat/Seurat.rds")
class(pbmc)
##将seurat格式转换成SingleCellExperiment格式。
pbmc.sce <- as.SingleCellExperiment(pbmc)
##这个分组就很随意，定义前4527个细胞是第一类，剩下的是第二类。
group <- factor(c(rep(1,4527),rep(2,9054)))
##计算差异表达基因即可，但是这个需要很久、很久、很久。。。。
results <- DEsingle(counts = pbmc.sce, group = group)
得到结果之后，就可以开始进行可视化展示了。
关于绘图的部分，我还需要再摸索下。
稍安勿躁，会持续更新。。。。。。
library(pheatmap)
library(ggplot2)
library(ggrepel)
library(export)
