單細(xì)胞轉(zhuǎn)錄調(diào)控網(wǎng)絡(luò)分析工具（TF，motifs）-SCENIC

相信很多專題應(yīng)該已經(jīng)講過SCENIC的使用，這兒，針對我自己對SCENIC的學(xué)習(xí)使用，遇到的問題做一個詳細(xì)的講解說明，希望能對大家理解學(xué)習(xí)應(yīng)用SCENIC起到一些幫助。

SCENIC (Single-Cell rEgulatory Network Inference and Clustering) 是從單細(xì)胞RNA數(shù)據(jù)推斷基因調(diào)控網(wǎng)絡(luò)及其相關(guān)細(xì)胞狀態(tài)的工具。
作者將SCENIC應(yīng)用于腫瘤和小鼠大腦單細(xì)胞圖譜數(shù)據(jù)，證明了順式調(diào)控網(wǎng)絡(luò)分析能夠有助于深入挖掘細(xì)胞異質(zhì)性背后的生物學(xué)意義，并為疾病的診斷、治療以及發(fā)育分化的研究提供有價值的線索。

SCENIC在2017年首先發(fā)表于nature methods，2020年又將流程整理后發(fā)表于nature protocls。有R版和python版本。
SCENIC: single-cell regulatory network inference and clustering https://www.nature.com/articles/nmeth.4463
A scalable SCENIC workflow for single-cell gene regulatory network analysis https://www.nature.com/articles/s41596-020-0336-2
這兒主要講述的是R版（無論哪種語言，代碼都可以copy到，重要的是理解做的意義）

作者提供了SCENIC 詳細(xì)的工作流程：http://htmlpreview.github.io/?https://github.com/aertslab/SCENIC/blob/master/inst/doc/SCENIC_Running.html
構(gòu)建基因調(diào)控網(wǎng)絡(luò)gene regulatory network (GRN)：
1），根據(jù)共表達(dá)確定每個 TF 的潛在靶基因：過濾表達(dá)式矩陣并運(yùn)行 GENIE3/GRNBoost；
基于共表達(dá)推斷轉(zhuǎn)錄因子與候選靶基因之間的共表達(dá)模塊。
2），使用RcisTarget對每個共表達(dá)模塊進(jìn)行順式調(diào)控基序（cis-regulatory motif）分析
3），使用AUCell算法對細(xì)胞中的regulon活性進(jìn)行評分（計算 AUC）
4），根據(jù)regulons活性確定穩(wěn)定的細(xì)胞狀態(tài)并探索結(jié)果
簡之，基于共表達(dá)和DNA?；?（motif）分析推斷基因調(diào)控網(wǎng)絡(luò) ，然后在每個細(xì)胞中分析網(wǎng)絡(luò)活性以鑒定細(xì)胞狀態(tài)。

安裝相關(guān)包
可參照官網(wǎng)：http://htmlpreview.github.io/?https://github.com/aertslab/SCENIC/blob/master/inst/doc/SCENIC_Setup.html

if (!requireNamespace("BiocManager", quietly = TRUE)) install.packages("BiocManager")
BiocManager::install(c("AUCell", "RcisTarget"))
BiocManager::install(c("GENIE3")) # Optional. Can be replaced by GRNBoost
BiocManager::install(c("zoo", "mixtools", "rbokeh"))
BiocManager::install(c("DT", "NMF", "ComplexHeatmap", "R2HTML", "Rtsne"))
BiocManager::install(c("doMC", "doRNG"))
if (!requireNamespace("devtools", quietly = TRUE)) install.packages("devtools")
devtools::install_github("aertslab/SCopeLoomR", build_vignettes = TRUE)

devtools::install_github("aertslab/SCENIC")  ##導(dǎo)入SCENIC
#這兒我把SCENIC安裝在服務(wù)器上，安裝過程中，也遇到了好多錯誤，主要是一些包導(dǎo)入錯誤（主要是網(wǎng)絡(luò)和環(huán)境引起的），把錯誤的包重新用conda進(jìn)行安裝，一定要耐心哦。

物種特異數(shù)據(jù)庫（人，鼠，果蠅）

#RcisTarget 的物種特定數(shù)據(jù)庫
dbFiles1 <- c("https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc9nr/gene_based/hg19-500bp-upstream-7species.mc9nr.feather",
"https://resources.aertslab.org/cistarget/databases/homo_sapiens/hg19/refseq_r45/mc9nr/gene_based/hg19-tss-centered-10kb-7species.mc9nr.feather")   
#人， SCENIC 使用數(shù)據(jù)庫對基因啟動子（TSS 上游高達(dá) 500bp）和 TSS 周圍 20kb (+/-10kbp) 中的motifs.
dbFiles2 <- c("https://resources.aertslab.org/cistarget/databases/mus_musculus/mm9/refseq_r45/mc9nr/gene_based/mm9-500bp-upstream-7species.mc9nr.feather",
"https://resources.aertslab.org/cistarget/databases/mus_musculus/mm9/refseq_r45/mc9nr/gene_based/mm9-tss-centered-10kb-7species.mc9nr.feather")  #鼠
dbFiles3 <-c("https://resources.aertslab.org/cistarget/databases/drosophila_melanogaster/dm6/flybase_r6.02/mc8nr/gene_based/dm6-5kb-upstream-full-tx-11species.mc8nr.feather") #果蠅
dir.create("./SCENIC") #構(gòu)建SCENIC分析文件夾
setwd("SCENIC")
dir.create("./SCENIC/cisTarget_databases")#構(gòu)建數(shù)據(jù)庫文件夾
setwd("cisTarget_databases")
for(featherURL in c(dbFiles1,dbFiles2,dbFiles3))
{
  download.file(featherURL, destfile=basename(featherURL)) # saved in current dir
}

導(dǎo)入包，進(jìn)行SEUNIC分析
這兒我們?yōu)榱私y(tǒng)一，也為了方便，還是使用之前簡書提到的pbmc3k數(shù)據(jù)集，畢竟目前使用10x數(shù)據(jù)的人比較多。當(dāng)然其它數(shù)據(jù)格式的導(dǎo)入作者也提供了代碼接入。作者也提供了一個例子數(shù)據(jù)集作為輸入，我在文章最后提供對應(yīng)代碼。
pbmc3k數(shù)據(jù)集相關(guān)下載，seurat聚類都可參照前面的簡書：http://m.itdecent.cn/p/adda4536b2cb

#導(dǎo)入相關(guān)R包
library(Seurat)
library(tidyverse)
library(patchwork)
library(SCENIC)

導(dǎo)入pbmc3k數(shù)據(jù)，并對細(xì)胞類型定義。（之前只是對pbmc3k數(shù)據(jù)進(jìn)行了聚類，但是并未定義細(xì)胞類型，這兒對細(xì)胞類型進(jìn)行定義）

pbmc <-readRDS("/home/wucheng/jianshu/function/data/pbmc.rds") 
setwd("/home/wucheng/SCENIC/SCENIC_pbmc3k") 
current.cluster.ids <- c(0:8)
new.cluster.ids <- c("Naive CD4 T", "CD14+ Mono", "Memory CD4 T", "B", "CD8 T", "FCGR3A+ Mono", "NK", "DC", "Platelet")
pbmc@meta.data$celltype <- plyr::mapvalues(x = pbmc@meta.data[,"seurat_clusters"], from = current.cluster.ids, to = new.cluster.ids)
head(pbmc@meta.data)
pdf(file="celltype_umap.pdf",width=10,height=10)    
DimPlot(pbmc, reduction = "umap", group.by= "celltype",label = TRUE, pt.size = 0.5) + NoLegend()
dev.off()
pbmc3k_celltype

準(zhǔn)備細(xì)胞meta信息

cellInfo <- data.frame(pbmc@meta.data)
colnames(cellInfo)[which(colnames(cellInfo)=="orig.ident")] <- "sample"
colnames(cellInfo)[which(colnames(cellInfo)=="nFeature_RNA")] <- "nGene"
colnames(cellInfo)[which(colnames(cellInfo)=="nCount_RNA")] <- "nUMI"
colnames(cellInfo)[which(colnames(cellInfo)=="seurat_clusters")] <- "cluster"
colnames(cellInfo)[which(colnames(cellInfo)=="celltype")] <- "celltype"
cellInfo <- cellInfo[,c("sample","nGene","nUMI","cluster","celltype")]
saveRDS(scenicOptions, file="int/cellInfo.Rds")
#設(shè)置celltype顏色，也可不設(shè)置，使用默認(rèn)顏色
colVars <- list(celltype=c("Naive CD4 T"="forestgreen",  "CD14+ Mono"="darkorange",  "Memory CD4 T"="magenta4", 
"B"="hotpink",  "CD8 T"="red3",  "FCGR3A+ Mono"="skyblue",   "NK"="darkblue",  "DC"="yellow",  "Platelet"="grey"))
colVars$celltype <- colVars$celltype[intersect(names(colVars$celltype), cellInfo$celltype)]
saveRDS(colVars, file="int/colVars.Rds")

初始化 SCENIC 設(shè)置,設(shè)置分析環(huán)境

org <- "hgnc" #mgi 代表鼠標(biāo)，hgnc 代表人，dmel 代表蒼蠅 ，我們這兒pbmc3k為人的，選擇hgnc
dbDir <- "./SCENIC/cisTarget_databases" # RcisTarget databases location
myDatasetTitle <- "SCENIC example on pbmc3k" # choose a name for your analysis
data(defaultDbNames)
dbs <- defaultDbNames[[org]]
scenicOptions <- initializeScenic(org=org, dbDir=dbDir, dbs=dbs, datasetTitle=myDatasetTitle, nCores=10)  #nCores=10，代表開啟10個線程計算
scenicOptions@inputDatasetInfo$cellInfo <- "int/cellInfo.Rds"
scenicOptions@inputDatasetInfo$colVars <- "int/colVars.Rds"
saveRDS(scenicOptions, file="int/scenicOptions.Rds")

構(gòu)建Co-expression network

exprMat <- as.matrix(pbmc@assays$RNA@counts) #準(zhǔn)備表達(dá)矩陣，為了節(jié)省計算資源，隨機(jī)只抽取部分細(xì)胞，計算共表達(dá)網(wǎng)絡(luò)
genesKept <- geneFiltering(exprMat, scenicOptions=scenicOptions,minCountsPerGene=3*.01*ncol(exprMat),minSamples=ncol(exprMat)*.01) #基因過濾/選擇，去除最有可能是噪音的基因
exprMat_filtered <- exprMat[genesKept, ]
dim(exprMat_filtered)
runCorrelation(exprMat_filtered, scenicOptions) ##計算相關(guān)性矩陣，1.2_corrMat.Rds：基因之間的相關(guān)性矩陣

根據(jù)表達(dá)數(shù)據(jù)推斷潛在的轉(zhuǎn)錄因子靶標(biāo)，使用 GENIE3 或 GRNBoost，GENIE3 非常耗時且計算量大（在 3-5k 單元的數(shù)據(jù)集上需要幾個小時或幾天的時間），GRNboost可在很短的時間內(nèi)提供與 GENIE3 類似的結(jié)果，這兒使用的R，選擇GENIC3

##GENIE3，GENIE3 的輸入通常是一個表達(dá)式矩陣和一個候選調(diào)節(jié)器列表。
exprMat_filtered <- log2(exprMat_filtered+1)  #完整矩陣
runGenie3(exprMat_filtered, scenicOptions, nParts = 10) #nParts參數(shù)，是把表達(dá)矩陣分成n份分開計算
#1.4，GENIE3_linkList.Rds：GENIE3最終結(jié)果，是把“1.3_”開頭的文件合并在一起

運(yùn)行到這兒，才進(jìn)入SCENIC關(guān)鍵分析步驟
1. 獲取共表達(dá)模塊
2. 獲取調(diào)節(jié)子（使用 RcisTarget）：TF 基序分析）
3. 對細(xì)胞中的 GRN評分（使用 AUCell）
4. 根據(jù) GRN 活動對細(xì)胞進(jìn)行聚類

exprMat_log <- log2(exprMat+1) #log標(biāo)準(zhǔn)話原始矩陣
#scenicOptions <- readRDS("int/scenicOptions.Rds")
scenicOptions@settings$verbose <- TRUE
scenicOptions@settings$nCores <- 1
scenicOptions@settings$seed <- 123
scenicOptions <- runSCENIC_1_coexNetwork2modules(scenicOptions) #1. 獲取共表達(dá)模塊
scenicOptions <- runSCENIC_2_createRegulons(scenicOptions)  #2. 獲取regulons
#可使用參數(shù)coexMethod=c("top5perTarget")，coexMethod可供選擇，作者嘗試了多種策略過濾低相關(guān)性TF-Target,并建議是6種過濾標(biāo)準(zhǔn)都用,其實(shí)也并沒有耗時多少。默認(rèn)都計算。
#w001：以每個TF為核心保留weight>0.001的基因形成共表達(dá)模塊；#w005，#top50
#top5perTarget：每個基因保留weight值top5的TF得到精簡的TF-Target關(guān)聯(lián)表，然后把基因分配給TF構(gòu)建共表達(dá)模塊；#top10perTarget；#top50perTarget
scenicOptions <- runSCENIC_3_scoreCells(scenicOptions, exprMat_log) #3. 對細(xì)胞中的 GRN（調(diào)節(jié)子）評分
saveRDS(scenicOptions, file="int/scenicOptions.Rds") 
runSCENIC_4_aucell_binarize(scenicOptions, exprMat=exprMat_log)

運(yùn)行完runSCENIC_2_createRegulons后產(chǎn)生：
output/Step2_MotifEnrichment.tsv
output/Step2_MotifEnrichment_preview.html
output/Step2_regulonTargetsInfo.tsv 文件

Step2_MotifEnrichment.tsv打開如下：


MotifEnrichment.tsv

geneSet:基因集的名稱；motif: motif ID,NES：基因集中motif的標(biāo)準(zhǔn)化富集分?jǐn)?shù)；
AUC：曲線下面積（用于計算 NES）； TFinDB：指示突出顯示的TF是否包含在高置信度注釋（兩個星號）或低置信度注釋（一個星號）中；TF_highConf：根據(jù)“motifAnnot_highConfCat”注釋到motif的轉(zhuǎn)錄因子；TF_lowConf：根據(jù)“motifAnnot_lowConfCat”注釋到motif的轉(zhuǎn)錄因子；richedGenes：在給定基序上排名較高的基因；nErnGenes：排名靠前的基因數(shù)量；rankAtMax：在最大富集處排序，用于確定富集基因的數(shù)量。
Step2_MotifEnrichment_preview.html


MotifEnrichment_preview.html

顯著富集的motif注釋信息，表頭與上個文件相同。
Step2_regulonTargetsInfo.tsv
Step2_regulonTargetsInfo.tsv

TF：轉(zhuǎn)錄因子名稱；gene：TF靶基因名稱；nMotif：靶基因在數(shù)據(jù)庫的motif數(shù)量；bestMotif：最顯著富集的motif名稱
NES：標(biāo)準(zhǔn)富集分?jǐn)?shù)，分值越高越顯著；highConfAnnot：是不是高可信注釋；Genie3Weight：TF與靶基因的相關(guān)性權(quán)重

運(yùn)行完runSCENIC_3_scoreCells后，產(chǎn)生regulonAUC矩陣，可用熱圖展示，t-SNE圖分別展示每個regulon的活性分布情況。


Step3_RegulonActivity_heatmap

runSCENIC_4_aucell_binarize，是進(jìn)行 二進(jìn)制轉(zhuǎn)換及衍生分析 ，這步可不做


Step4_BinaryRegulonActivity_Heatmap_all

這兒，也可使用GRNs在每個細(xì)胞中的活性得分，來對細(xì)胞進(jìn)行聚類，或許會得到不一樣的細(xì)胞狀態(tài)分類

nPcs <- c(5,15,50)
scenicOptions@settings$seed <- 123 # same seed for all of them
# Run t-SNE with different settings:
fileNames <- tsneAUC(scenicOptions, aucType="AUC", nPcs=nPcs, perpl=c(5,15,50), onlyHighConf=TRUE, filePrefix="int/tSNE_oHC")
# Plot as pdf (individual files in int/):
pdf("int/AUC_tsne.pdf",width=10,height=10)
plotTsne_compareSettings(fileNames, scenicOptions, showLegend=FALSE, varName="celltype", cex=.5)
dev.off()
AUC_tsne

當(dāng)然，其實(shí)分析到這兒還是有點(diǎn)懵，regulons可能需要進(jìn)一步挑選，可視化，結(jié)合Seurat聚類結(jié)果分析比較簡單且有實(shí)際意義

pbmc <-readRDS("/home/wucheng/jianshu/function/data/pbmc.rds")  #導(dǎo)入pbmc3k數(shù)據(jù)
current.cluster.ids <- c(0:8)
new.cluster.ids <- c("Naive CD4 T", "CD14+ Mono", "Memory CD4 T", "B", "CD8 T", "FCGR3A+ Mono", "NK", "DC", "Platelet")
pbmc@meta.data$celltype <- plyr::mapvalues(x = pbmc@meta.data[,"seurat_clusters"], from = current.cluster.ids, to = new.cluster.ids)
head(pbmc@meta.data)

##導(dǎo)入原始regulonAUC矩陣，也就是運(yùn)行完runSCENIC_3_scoreCells后產(chǎn)生的AUC矩陣，可以查看每個GRN在每個細(xì)胞中的AUC活性打分
setwd("/home/wucheng/SCENIC/SCENIC_pbmc3k")
AUCmatrix <- readRDS("int/3.4_regulonAUC.Rds")
AUCmatrix <- data.frame(t(AUCmatrix@assays@data@listData$AUC), check.names=F)
RegulonName_AUC <- colnames(AUCmatrix)
RegulonName_AUC <- gsub(' \\(','_',RegulonName_AUC) #把(替換成_
RegulonName_AUC <- gsub('\\)','',RegulonName_AUC) #把)去掉，最后如把 KLF3_extended (79g) 替換成KLF3_extended_79g
colnames(AUCmatrix) <- RegulonName_AUC
pbmcauc <- AddMetaData(pbmc, AUCmatrix) #把AUC矩陣添加到pbmc的metadata信息中
pbmcauc@assays$integrated <- NULL
saveRDS(pbmcauc,'pbmcauc.rds')

##導(dǎo)入二進(jìn)制regulonAUC矩陣
BINmatrix <- readRDS("int/4.1_binaryRegulonActivity.Rds")
BINmatrix <- data.frame(t(BINmatrix), check.names=F)
RegulonName_BIN <- colnames(BINmatrix)
RegulonName_BIN <- gsub(' \\(','_',RegulonName_BIN)
RegulonName_BIN <- gsub('\\)','',RegulonName_BIN)
colnames(BINmatrix) <- RegulonName_BIN
pbmcbin <- AddMetaData(pbmc, BINmatrix)
pbmcbin@assays$integrated <- NULL
saveRDS(pbmcbin, 'pbmcbin.rds')

##利用Seurat可視化AUC
dir.create('scenic_seurat')
#FeaturePlot
GRNs <-intersect(colnames(AUCmatrix),colnames(BINmatrix))
for(i in 1:length(GRNs)){
p1 = FeaturePlot(pbmcauc, features=GRNs[i], label=T, reduction = 'umap')
p2 = FeaturePlot(pbmcbin, features=GRNs[i], label=T, reduction = 'umap')
p3 = DimPlot(pbmc, reduction = 'umap', group.by = "celltype", label=T)
plotc = p1|p2|p3
ggsave(paste("scenic_seurat/",GRNs[i],".png",sep=""), plotc, width=14 ,height=4)
}
YY1_512g

CEBPA_31g

可以發(fā)現(xiàn)，有一些轉(zhuǎn)錄因子行成的調(diào)控網(wǎng)絡(luò)主要調(diào)控某些特異細(xì)胞類型，比如這兒的CEBPA。

#RidgePlot&VlnPlot ,和單細(xì)胞展示基因表達(dá)一樣，這兒可以利用小提琴圖展示GRNs的活性得分
for(i in 1:length(GRNs)){
p1 = RidgePlot(pbmcauc, features =GRNs[i], group.by="celltype") + theme(legend.position='none')
p2 = VlnPlot(pbmcauc, features =GRNs[i], pt.size = 0, group.by="celltype") + theme(legend.position='none')
plotc = p1 + p2
ggsave(paste("scenic_seurat/","Ridge-Vln_",GRNs[i],".png",sep=""),plotc, width=10, height=8)
}
CEBPA_31g

當(dāng)然也可以用熱圖可視化SCENIC結(jié)果

library(pheatmap)
cellInfo <- readRDS("int/cellInfo.Rds")
celltype = as.data.frame(subset(cellInfo,select = 'celltype'))
AUCmatrix <- t(AUCmatrix)
BINmatrix <- t(BINmatrix)
#挑選部分感興趣的regulons
my.regulons <- intersect(rownames(AUCmatrix),rownames(BINmatrix))
myAUCmatrix <- t(AUCmatrix[rownames(AUCmatrix)%in%my.regulons,])
myBINmatrix <- t(BINmatrix[rownames(BINmatrix)%in%my.regulons,])
#使用regulon原始AUC值繪制熱圖
pdf("scenic_seurat/AUC_heatmap.pdf",width=10,height=20)
pheatmap(myAUCmatrix, show_colnames=F, annotation_col=celltype)
dev.off()
#使用regulon二進(jìn)制AUC值繪制熱圖
pdf("scenic_seurat/BIN_heatmap.pdf",width=10,height=20)
pheatmap(myBINmatrix, show_colnames=F, annotation_col=celltype,color = colorRampPalette(colors = c("white","black"))(100))
dev.off()
AUC_heatmap-BIN_heatmap

如果覺得regulons太多了，有一些是沒有意義的，可以使用RSS來識別挑選各個細(xì)胞類型代表性regulons，進(jìn)行分析

regulonAUC <-loadInt(scenicOptions, "aucell_regulonAUC")
rss <- calcRSS(AUC=getAUC(regulonAUC), cellAnnotation=cellInfo[colnames(regulonAUC), "celltype"])
rssPlot <-plotRSS(rss) #大小 rss評分，顏色 Z-score 
ggsave('int/rss.png', rssPlot$plot, width=5 ,height=20)
rss

如果不想用pbmc3k數(shù)據(jù)，也可使用作者示例的小鼠大腦單細(xì)胞數(shù)據(jù)進(jìn)行SCENIC分析

library(Seurat)
library(tidyverse)
library(patchwork)
library(SCENIC)
library(SCopeLoomR)
dir.create("./SCENIC_MouseBrain") 
setwd("SCENIC_MouseBrain")
loomPath <- system.file(package="SCENIC", "examples/mouseBrain_toy.loom") ##作者提供的實(shí)例數(shù)據(jù)
loom <- open_loom(loomPath)
exprMat <- get_dgem(loom)
cellInfo <- get_cell_annotation(loom)
close_loom(loom)
dim(exprMat)
cellInfo <- data.frame(cellInfo)
cbind(table(cellInfo$CellType))
##初始化 SCENIC 設(shè)置
org <- "mgi" #mgi 代表鼠標(biāo)，hgnc 代表人，dmel 代表蒼蠅
dbDir <- "/home/wucheng/SCENIC/cisTarget_databases" # RcisTarget databases location
myDatasetTitle <- "SCENIC example on Mouse brain" # choose a name for your analysis
data(defaultDbNames)
dbs <- defaultDbNames[[org]]
scenicOptions <- initializeScenic(org=org, dbDir=dbDir, dbs=dbs, datasetTitle=myDatasetTitle, nCores=10) 
scenicOptions@inputDatasetInfo$cellInfo <- "int/cellInfo.Rds"
scenicOptions@inputDatasetInfo$colVars <- "int/colVars.Rds"
saveRDS(scenicOptions, file="int/scenicOptions.Rds") 
##構(gòu)建Co-expression network
#1，根據(jù)表達(dá)數(shù)據(jù)推斷潛在的轉(zhuǎn)錄因子靶標(biāo)，使用 GENIE3 或 GRNBoost
#GENIE3 非常耗時且計算量大（在 3-5k 單元的數(shù)據(jù)集上需要幾個小時或幾天的時間）
#GRNboost可在很短的時間內(nèi)提供與 GENIE3 類似的結(jié)果
#基因過濾/選擇，去除最有可能是噪音的基因
genesKept <- geneFiltering(exprMat, scenicOptions=scenicOptions,minCountsPerGene=3*.01*ncol(exprMat)，minSamples=ncol(exprMat)*.01)
runCorrelation(exprMat_filtered, scenicOptions) #GENIE3，GENIE3 的輸入通常是一個表達(dá)式矩陣和一個候選調(diào)節(jié)器列表。
exprMat_filtered <- log2(exprMat_filtered+1) 
runGenie3(exprMat_filtered, scenicOptions, nParts = 10)
#構(gòu)建 GRN 并為其評分#1. 獲取共表達(dá)模塊 
2. 獲取調(diào)節(jié)子（使用 RcisTarget）：TF 基序分析）
#識別細(xì)胞狀態(tài)： 3. 對細(xì)胞中的 GRN（調(diào)節(jié)子）評分（使用 AUCell） 
#4. 根據(jù) GRN 活動對細(xì)胞進(jìn)行聚類
exprMat_log <- log2(exprMat+1)
library(SCENIC)
scenicOptions <- readRDS("int/scenicOptions.Rds")
scenicOptions@settings$verbose <- TRUE
scenicOptions@settings$nCores <- 1
scenicOptions@settings$seed <- 123
scenicOptions <- runSCENIC_1_coexNetwork2modules(scenicOptions)
scenicOptions <- runSCENIC_2_createRegulons(scenicOptions, coexMethod=c("top5perTarget")) 
scenicOptions <- runSCENIC_3_scoreCells(scenicOptions, exprMat_log)
saveRDS(scenicOptions, file="int/scenicOptions.Rds")
AUC

也可參照其他優(yōu)質(zhì)簡書文章對SECNIC的理解，如下：
https://cloud.tencent.com/developer/article/1692240
http://m.itdecent.cn/p/bd0df23fbd60
https://www.cnblogs.com/raisok/p/12425225.html

希望大家關(guān)注點(diǎn)贊，謝謝?。。。。。。。。。。?！
有不清楚的可以一起探討，謝謝?。。?！

?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成，瀏覽時請結(jié)合常識與多方信息審慎甄別。
平臺聲明：文章內(nèi)容（如有圖片或視頻亦包括在內(nèi)）由作者上傳并發(fā)布，文章內(nèi)容僅代表作者本人觀點(diǎn)，簡書系信息發(fā)布平臺，僅提供信息存儲服務(wù)。

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容