Skip to content

Commit 3416ef1

Browse files
WWXkenmoWWXkenmo
authored
upload NetID v0.1.1
1 parent 6ad1ffd commit 3416ef1

File tree

4 files changed

+498
-433
lines changed

4 files changed

+498
-433
lines changed

‎DESCRIPTION

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -2,11 +2,11 @@ Package: NetID
22
Type: Package
33
Title: A Meta-cell Base Method to Infer Fate Specific Gene Regulatory
44
Network from Single Cell Gene Expression Data
5-
Version: 0.1.0
5+
Version: 0.1.1
66
biocViews:
77
Imports: Matrix, lmtest, glmnet, RobustRankAggreg, irlba, rsvd,
88
reticulate, Seurat, Hmisc, mclust, pracma, doParallel, doRNG,
9-
RaceID (>= 0.2.9), igraph, rARPACK
9+
RaceID (== 0.3.9), igraph, rARPACK
1010
Suggests: GENIE3, SingleCellExperiment, knitr, rmarkdown
1111
VignetteBuilder: knitr
1212
Author: Weixu Wang and Dominic Grün

‎R/FateDynamic.R

Lines changed: 147 additions & 110 deletions
Original file line numberDiff line numberDiff line change
@@ -21,121 +21,158 @@
2121
#'
2222
#' @export
2323
#'
24-
FateDynamic <- function(sce,global_params = list(),palantir_params = list(),cellrank_params = list(),work_dir = NULL){
25-
## Build Annodata object using SingleCellExperiments
26-
if(!is.null(work_dir)){setwd(work_dir)}
27-
env = environment()
28-
29-
anndata <- reticulate::import('anndata', convert = FALSE)
30-
sc <- reticulate::import('scanpy', convert = FALSE)
31-
scv <- reticulate::import('scvelo', convert = FALSE)
32-
33-
sce <- sce[,colSums(SummarizedExperiment::assays(sce)$spliced)>0]
34-
X <- SummarizedExperiment::assays(sce)$spliced
35-
obs <- as.data.frame(SummarizedExperiment::colData(sce))
36-
var <- data.frame(genes = rownames(sce))
37-
X <- X[!duplicated(var$genes),]
38-
var <- data.frame(genes = rownames(X))
39-
rownames(var) <- rownames(X)
40-
obsm <- NULL
41-
reductions <- SingleCellExperiment::reducedDimNames(sce)
42-
if (length(reductions) > 0) {
43-
obsm <- sapply(
44-
reductions,
45-
function(name) as.matrix(SingleCellExperiment::reducedDim(sce, name)),
46-
simplify = FALSE
47-
)
48-
names(obsm) <- paste0('X_', tolower(SingleCellExperiment::reducedDimNames(sce)))
49-
}
50-
51-
layers <- list()
52-
for (layer in c("spliced","unspliced")) {
53-
mat <- SummarizedExperiment::assays(sce)[[layer]]
54-
mat <- mat[rownames(X),]
55-
layers[[layer]] <- Matrix::t(mat)
56-
}
57-
58-
adata <- anndata$AnnData(
59-
X = Matrix::t(X),
60-
obs = obs,
61-
var = var,
62-
obsm = obsm,
63-
layers = layers
64-
)
65-
66-
adata$write("adata.h5ad")
67-
####
68-
# input global parameters
69-
global_params <- check_global_params(global_params)
70-
list2env(global_params,env)
71-
## check parameters
72-
palantir_params <- check_palantir_params(palantir_params)
73-
list2env(palantir_params,env)
74-
cellrank_params <- check_cellrank_params(cellrank_params)
75-
list2env(cellrank_params,env)
24+
FateDynamic <- function (sce, global_params = list(), palantir_params = list(),
25+
cellrank_params = list(), work_dir = NULL)
26+
{
27+
if (!is.null(work_dir)) {
28+
setwd(work_dir)
29+
}
30+
env = environment()
31+
anndata <- reticulate::import("anndata", convert = FALSE)
32+
sc <- reticulate::import("scanpy", convert = FALSE)
33+
scv <- reticulate::import("scvelo", convert = FALSE)
7634

77-
# input parameters
78-
writeLines(paste("Using palantir to perform cell fate analysis...",sep=""))
79-
n.cores = 8
80-
if(is.null(palantir_params[["start_cell"]])){
81-
if(method == "cytotrace"){
82-
writeLines("Using CyToTRACE to identify start cell...")
83-
if (require(CytoTRACE)) {
84-
res = CytoTRACE::CytoTRACE(as.matrix(X),ncores = n.cores)$CytoTRACE
85-
res = as.matrix(res)
86-
}
35+
if ( class(sce) == "SingleCellExperiment" | class(sce) == "Seurat" ){
36+
if ( class(sce) == "Seurat" ) sce <- as.SingleCellExperiment(sce)
37+
if ( ! ( "spliced" %in% names(assays(sce)) & "unspliced" %in% names(assays(sce)) ) ){
38+
assaySp = "counts"
39+
assayUn = NA
40+
}else{
41+
assaySp = "spliced"
42+
assayUn = "unspliced"
43+
}
44+
45+
46+
sce <- sce[, colSums(SummarizedExperiment::assays(sce)[[assaySp]]) > 0]
47+
X <- SummarizedExperiment::assays(sce)[[assaySp]]
48+
obs <- as.data.frame(SummarizedExperiment::colData(sce))
49+
var <- data.frame(genes = rownames(sce))
50+
X <- X[!duplicated(var$genes), ]
51+
var <- data.frame(genes = rownames(X))
52+
rownames(var) <- rownames(X)
53+
obsm <- NULL
54+
reductions <- SingleCellExperiment::reducedDimNames(sce)
55+
56+
if (length(reductions) > 0) {
57+
obsm <- sapply(reductions, function(name) as.matrix(SingleCellExperiment::reducedDim(sce,name)), simplify = FALSE)
58+
names(obsm) <- paste0("X_", tolower(SingleCellExperiment::reducedDimNames(sce)))
59+
}
60+
layers <- list()
61+
for (layer in c(assaySp, assayUn)) {
62+
if ( !is.na(layer) ){
63+
mat <- SummarizedExperiment::assays(sce)[[layer]]
64+
mat <- mat[rownames(X), ]
65+
layers[[layer]] <- Matrix::t(mat)
66+
}
67+
}
68+
}else if (class(sce) == "SCseq" ) {
69+
assaySp = "counts"
70+
assayUn = NA
71+
72+
X <- getExpData(sce,genes=rownames(sce@expdata))
73+
74+
obs <- data.frame(cpart= as.character(sce@cpart))
75+
rownames(obs) <- names(sce@cpart)
76+
var <- data.frame(genes = rownames(X))
77+
rownames(var) <- rownames(X)
78+
obsm <- NULL
79+
layers <- list()
80+
layers[[assaySp]] <- Matrix::t(X)
8781
}
88-
if(method == "scent"){
89-
writeLines("Using SCENT to identify start cell...")
90-
object <- Seurat::CreateSeuratObject(X)
91-
object <- Seurat::NormalizeData(object)
92-
X_norm <- Seurat::GetAssayData(object)
93-
rm(object);gc()
94-
if(length(intersect(list.files(getwd()),"PPI.rds")) == 1){
95-
writeLines("Find PPI object at local dictionary, Read PPI...")
96-
PPI <- readRDS("PPI.rds")
97-
}else{
98-
PPI <- getPPI_String(X_norm,species=species)
99-
}
100-
if (require(SCENT)) {
101-
res = SCENT::CompCCAT(X_norm, PPI)
102-
}
103-
saveRDS(PPI,file="PPI.rds")
104-
res = as.matrix(res)
105-
rownames(res) <- colnames(X_norm)
106-
rm(X_norm,PPI);gc()
82+
adata <- anndata$AnnData(X = Matrix::t(X), obs = obs, var = var, obsm = obsm, layers = layers)
83+
adata$write("adata.h5ad")
84+
global_params <- check_global_params(global_params)
85+
list2env(global_params, env)
86+
palantir_params <- check_palantir_params(palantir_params)
87+
list2env(palantir_params, env)
88+
cellrank_params <- check_cellrank_params(cellrank_params)
89+
list2env(cellrank_params, env)
90+
writeLines(paste("Using palantir to perform cell fate analysis...", sep = ""))
91+
n.cores = 8
92+
if (is.null(palantir_params[["start_cell"]])) {
93+
if (method == "cytotrace") {
94+
writeLines("Using CyToTRACE to identify start cell...")
95+
if (require(CytoTRACE)) {
96+
res = CytoTRACE::CytoTRACE(as.matrix(X), ncores = n.cores)$CytoTRACE
97+
res = as.matrix(res)
98+
}
99+
}
100+
if (method == "scent") {
101+
writeLines("Using SCENT to identify start cell...")
102+
object <- Seurat::CreateSeuratObject(X)
103+
object <- Seurat::NormalizeData(object)
104+
X_norm <- Seurat::GetAssayData(object)
105+
rm(object)
106+
gc()
107+
if (length(intersect(list.files(getwd()), "PPI.rds")) ==
108+
1) {
109+
writeLines("Find PPI object at local dictionary, Read PPI...")
110+
PPI <- readRDS("PPI.rds")
111+
} else {
112+
PPI <- getPPI_String(X_norm, species = species)
113+
}
114+
if (require(SCENT)) {
115+
res = SCENT::CompCCAT(X_norm, PPI)
116+
}
117+
saveRDS(PPI, file = "PPI.rds")
118+
res = as.matrix(res)
119+
rownames(res) <- colnames(X_norm)
120+
rm(X_norm, PPI)
121+
gc()
122+
}
123+
if (method == "markergene") {
124+
res = as.matrix(X[root_gene, ])
125+
}
126+
start_cell = rownames(res)[which.max(res)]
107127
}
108-
if(method == "markergene"){
109-
res = as.matrix(X[root_gene,])
128+
terminalstate = terminal_state[1]
129+
if (is.null(terminalstate)) {
130+
terminalstate = "None"
131+
} else {
132+
for (i in 2:length(terminal_state)) {
133+
terminalstate = paste(terminalstate, " ", terminal_state[i],
134+
sep = "")
135+
}
110136
}
111-
start_cell = rownames(res)[which.max(res)]
112-
}
113-
## terminal_state
114-
terminalstate = terminal_state[1]
115-
if(is.null(terminalstate)){
116-
terminalstate = "None"
117-
}else{
118-
for(i in 2:length(terminal_state)){
119-
terminalstate = paste(terminalstate," ",terminal_state[i],sep="")
120-
}
121-
}
122-
123-
input = "adata.h5ad"
124-
fate_predict_py = Hmisc::capitalize(tolower(as.character(fate_predict)))
125-
plot = Hmisc::capitalize(tolower(as.character(plot)))
126-
#if(is.null(cluster_label)) cluster_label = ""
127-
script = system.file("/python/FateDynamic_py.py", package = "NetID")
128-
commandLines = paste('python ',script,' -i ',input,' -m ',min_counts,' -n ',nhvgs,' -pc ',npcs,' -k ', knn, ' -pl ',"True",' -dc ',ndcs,' -r ',start_cell,' -ts ',terminalstate,' -nw ',nwps,' -mo ',mode,' -p ',plot,' -f ',fate_predict_py, ' -c ', cluster_label, ' -w ',weight_connectivities,' -nc ', ncores, ' -t ',tolerance,sep="")
129-
#reticulate::py_run_file(system.file(commandLines, package = "NetID"))
130-
system(commandLines, wait=TRUE)
131-
132-
if(velo){
133-
writeLines(paste("Using cellrank with ",mode," mode velocity to perform cell fate analysis...",sep=""))
137+
input = "adata.h5ad"
138+
fate_predict_py = Hmisc::capitalize(tolower(as.character(fate_predict)))
139+
plot = Hmisc::capitalize(tolower(as.character(plot)))
134140
script = system.file("/python/FateDynamic_py.py", package = "NetID")
135-
commandLines = paste('python ',script,' -i ',input,' -m ',min_counts,' -n ',nhvgs,' -pc ',npcs,' -k ', knn, ' -pl ',"False",' -dc ',ndcs,' -r ',start_cell, ' -ts ',terminalstate, ' -nw ',nwps,' -mo ',mode,' -p ',plot,' -f ',fate_predict_py, ' -c ', cluster_label, ' -w ',weight_connectivities,' -nc ', ncores, ' -t ',tolerance,sep="")
136-
#reticulate::py_run_file(system.file(commandLines, package = "NetID"))
137-
system(commandLines, wait=TRUE)
138-
}
141+
if ( class(sce) == "SCseq" ){
142+
commandLines = paste("python ", script, " -i ", input, " -m ",
143+
min_counts, " -n ", nhvgs, " -pc ", npcs, " -k ", knn,
144+
" -pl ", "True", " -dc ", ndcs, " -r ", start_cell, " -ts ",
145+
terminalstate, " -nw ", nwps, " -mo ", mode, " -p ",
146+
plot, " -f ", fate_predict_py, " -c cpart",
147+
" -w ", weight_connectivities, " -nc ", ncores, " -t ",
148+
tolerance, sep = "")
149+
}else{
150+
commandLines = paste("python ", script, " -i ", input, " -m ",
151+
min_counts, " -n ", nhvgs, " -pc ", npcs, " -k ", knn,
152+
" -pl ", "True", " -dc ", ndcs, " -r ", start_cell, " -ts ",
153+
terminalstate, " -nw ", nwps, " -mo ", mode, " -p ",
154+
plot, " -f ", fate_predict_py, " -c ", cluster_label,
155+
" -w ", weight_connectivities, " -nc ", ncores, " -t ",
156+
tolerance, sep = "")
157+
}
158+
system(commandLines, wait = TRUE)
159+
160+
if (velo) {
161+
if (!is.na(assayUn) ){
162+
writeLines(paste("Using cellrank with ", mode, " mode velocity to perform cell fate analysis...", sep = ""))
163+
script = system.file("/python/FateDynamic_py.py", package = "NetID")
164+
commandLines = paste("python ", script, " -i ", input,
165+
" -m ", min_counts, " -n ", nhvgs, " -pc ", npcs,
166+
" -k ", knn, " -pl ", "False", " -dc ", ndcs, " -r ",
167+
start_cell, " -ts ", terminalstate, " -nw ", nwps,
168+
" -mo ", mode, " -p ", plot, " -f ", fate_predict_py,
169+
" -c ", cluster_label, " -w ", weight_connectivities,
170+
" -nc ", ncores, " -t ", tolerance, sep = "")
171+
system(commandLines, wait = TRUE)
172+
}else{
173+
writeLines(paste("No unspliced counts available! Falling back to Palantir...", sep = ""))
174+
}
175+
}
139176
}
140177

141178
getPPI_String <- function (object = NULL, species = 9606, score_threshold = 600,

0 commit comments

Comments
 (0)