单细胞代码解析-妇科癌症单细胞转录组及染色质可及性分析1：https://cloud.tencent.com/developer/article/2055573
单细胞代码解析-妇科癌症单细胞转录组及染色质可及性分析2：https://cloud.tencent.com/developer/article/2072069
这部分的是作者的PART2，我发现没有相关的注释文件基本上做不出来后面的东西，我只能一点一点的尝试读懂代码了。。。
这部分主要的是对两个去除双细胞的R包的代码进行解析。
图片.png
图片.png代码解析
代码语言：javascript复制###########################################################
# Part 2: Doublet detection and removal
###########################################################
#' # Doublet detection: 1) DoubletDecon 2) DoubletFinder 3) Take intersection of calls
######################################################################################
# # Doublet Decon
# # Add if else statement to regress out nCount RNA if needed before running DoubletDecon
##DoubletDecon：该方法结合反卷积分析（deconvolution）和鉴定细胞独有的基因表达来检测Doublet
#软件链接：[https://github.com/EDePasquale/DoubletDecon](https://github.com/EDePasquale/DoubletDecon)
library(DoubletDecon)
#length() 函数用于获取或设置向量(列表)或其他对象的长度；levels专门是处理factor变量的level属性用的；as.numeric：将因子变量（factor）转化为数值变量
idents.length <- length(levels(Idents(rna)))
for (i in levels(Idents(rna))){
i <- as.numeric(i)
levels(Idents(rna))[i] <- i -1
}
#Improved_Seurat_Pre_Process()
#as.factor () R语言中的函数用于将传递的对象 (通常是Vector)转换为Factor。
#Idents(rna) <- as.factor(Idents(rna))
seuratObject=rna
#Seurat创建对象和细胞过滤
newFiles=Improved_Seurat_Pre_Process(seuratObject, num_genes=50, write_files=F)
#dev.off()
#rawDataFile：表达矩阵（gene by cell）；groupsFile：包含组分配的文件名（3列：cell，组，组（数字或字符））；filename：唯一的文件名，输入文件的名字；location：应在其中存储输出的目录
#fullDataFile：包含完整表达式数据的文件名（gene by cell）。默认值为NULL。removeCC：通过KEGG富集去除细胞周期基因。默认值为FALSE。
#species：科学物种名称，KEGG ID，三个字母的物种缩写或NCBI ID。默认值为“ mmu”。rhop：平均值x中的x * SD以确定黑名单中相关性的上限。默认值为1。
#write：将输出文件写为.txt文件。默认值为TRUE。recluster：recluster反卷积使用Hopach或反卷积分类分别对doublet和非doublet进行分类。
#PMF：在双重确定标准中使用步骤3（独特的基因表达）。默认值为TRUE。useFull：使用完整的基因列表进行PMF分析。需要fullDataFile。默认值为FALSE。
#heatmap：是否生成热图的布尔值。默认值为TRUE。大于约3000个像元的数据集可能比较慢。重心：在解卷积中，将重心用作参考，而不是默认重心。
#num_doubs：用户定义的每对集群要生成的双峰数目。默认值为100。only50：仅使用由50％/ 50％的父单元格混合创建的合成对偶，而不是30％/ 70％和70％/ 30％的扩展选项，默认为FALSE。
#min_uniq：挽救群集所需的最小独特基因数，默认值为4。
results=Main_Doublet_Decon(rawDataFile=newFiles$newExpressionFile,
groupsFile=newFiles$newGroupsFile,
filename=SAMPLE.ID,
location="./DoubletDecon",
fullDataFile=NULL,
removeCC=FALSE,
species="hsa",
rhop=1.1,
write=F,
PMF=TRUE,
useFull=FALSE,
heatmap=FALSE,
centroids=TRUE,
num_doubs=100,
only50=FALSE,
min_uniq=4,
nCores=4)
#rownames获取和设置数据框架的行名；gsub ()函数是2R语言中处理正则表达式中的一种
decon.doublets <- rownames(results$Final_doublets_groups)
decon.doublets <- gsub("\.","-",decon.doublets)
# DoubletFinder
# Add modfieid Parameter sweep function to regress out nCOunt RNA if needed
# Add if else statement to regress out nCount RNA if needed before running DoubletFinder
#########################################################################################
#DoubletFinder是一种R包，可预测单细胞RNA 测序 数据中的doublet，具体解析[https://www.jianshu.com/p/b1947c4156ad](https://www.jianshu.com/p/b1947c4156ad)
#软件链接：[https://github.com/chris-mcginnis-ucsf/DoubletFinder](https://github.com/chris-mcginnis-ucsf/DoubletFinder)
library(DoubletFinder)
# Estimate Doublets
#################################################################################################
## pK Identification (no ground-truth) ---------------------------------------------------------------------------------------
#DoubletFinder 去除双细胞
sweep.res.list <- paramSweep_v3(rna, PCs = 1:50, sct = FALSE)
sweep.stats <- summarizeSweep(sweep.res.list , GT = FALSE)
#找最佳 PK
bcmvn<- find.pK(sweep.stats)
pK.1 <- as.numeric(unfactor(dplyr::arrange(bcmvn,desc(BCmetric))$pK[1]))
nExp_poi <- round(doublet.rate*length(colnames(counts))) ## Assuming 4.6% doublet formation rate - tailor for your dataset
# Run doubletFinder with pk.1
rna <- doubletFinder_v3(rna, PCs = 1:50, pN = 0.25, pK = pK.1, nExp = nExp_poi, reuse.pANN = FALSE, sct = FALSE)
meta <- rna@meta.data
#R 语言中的 paste0 () 函数用于连接所有没有分隔符的元素
doublet.column <- paste0("DF.classifications_0.25_",pK.1,"_",nExp_poi)
doublet.calls <- rna[[doublet.column]]
colnames(doublet.calls) <- "Call"
##dplyr::filter：数据过滤
rna.dub <- dplyr::filter(doublet.calls, Call == "Doublet")
rna.singlet <- dplyr::filter(doublet.calls, Call == "Singlet")
DF.doublets <- rownames(rna.dub)
# Intersect doublet calls
###########################################################################
head(DF.doublets)
head(decon.doublets)
#intersect () R语言中的函数用于查找两个对象的交集
doublets <- intersect(decon.doublets,DF.doublets)
rna@meta.data$Doublet.Call <- ifelse(rownames(rna@meta.data) %in% doublets,"TRUE","FALSE")
#FeatureScatte：在一组单个单元格中创建两个特征（通常是特征表达式）的散点图。细胞的颜色是由它们的身份类别决定的。皮尔逊两个特征之间的相关性显示在绘图上方。
FeatureScatter(rna,feature1 = "nCount_RNA",feature2 = "nFeature_RNA",group.by = "Doublet.Call") ggsave("Doublet_calls_FeatureScatter.png")
DimPlot(rna,group.by = "Doublet.Call") ggsave("Doublet_calls_UMAP.png")
saveRDS(doublets,"Intersection_DoubletDecon_DoubletFinder_doublets.rds")
cells <- colnames(rna)总结
作者这部分的内容主要是去了去除双细胞进行的分析，选用了两个R包DoubletDecon和DoubletFinder，然后对两个去除双细胞的结果进行相关性分析，去判断结果的可靠性。这里分析的主要目的是保证后面的分析是单细胞的结果，不会影响的后面的结果，因此在前面做了两重的分析，进行验证。
