Awesome
scLM
A R-based tool to do the automatic identification of co-expressed genes across mulitple single cell RNA-seq datasets simultaneously
1. Installation
devtools::install_github("QSong-github/scLM")
2. Input scRNA-seq Data File Format
scLM works with multiple single cell RNA-seq dataset as inputs. It also works with one single cell dataset. Bascially, the format looks like the following. Example data files can be found in the Data
folder.
CellID | Cell1ID | Cell2ID | Cell3ID | Cell4ID | ... |
---|---|---|---|---|---|
Gene1 | 12 | 0 | 0 | 0 | ... |
Gene2 | 125 | 0 | 298 | 0 | ... |
Gene3 | 0 | 0 | 0 | 0 | ... |
... | ... | ... | ... | ... | ... |
The gourd truth labels for cells in each dataset can also be input. The format is as following
Cell1ID | Lable1 |
---|---|
Cell2ID | Lable2 |
Cell3ID | Lable3 |
Cell4ID | Lable4 |
... | .. |
3. How to use
3.1 load scLM package
library(scLM)
recommend the linux system to run the codes
3.2 read in scRNA-seq datasets
#' load the example data
data("example1")
# or
# data("example2")
In this example, we define 3 co-expression clusters for each dataset of the input list example1 and example2
3.3 read in ground truth cell labels (this is optional)
data(example1.member)
# data(example2.member)
3.4 run the function with specific lambda across multiple datasets
result1 <- Multi_NB(datalist=example1, K=3, N=nrow(example1[[1]]))
result1 contains the latent variables accompanied with other coefficients, the identified co-expression clusters
3.5 evaluate of clustering results using ground truth (this is optional)
# calculate the Adjusted Rand Index
library(clues)
adjustedRand(result1$clusters,example1.member)
4. Examples and reproducible results
can be found using the example.R script
5. without prior knowledge, run the function with several lambda across multiple datasets
5.1 How to use
for (lambda in 1:20)
{
results <- Multi_NB(datalist=example1, K=lambda, N=nrow(example1[[1]]))
save(results,file=paste0('path1',lambda,'results.RData'))
}
Output results in a designated path "path1"
5.2 Identify the optimal lambda and co-expression clusters
files <- list.files(path=path1, pattern='results.RData')
# load all the results in the resA with the list structure
count <- 0
resA <- vector('list')
for ( i in files){
load(i)
count <- count +1
resA[[count]] <- results
names(resA[[count]]) <- paste0('res_',strsplit(strsplit(i,'_')[[1]][4],'results')[[1]][1])}
# identify the result with least BIC value
bicS <- lapply(1:length(resA),function(i){ Res <- resA[[i]][[1]]$BIC })
optimal.lambda <- grep(min(unlist(bicS)),unlist(bicS))
# optimal co-expression clusters
load(files[optimal.lambda])
opitmal.cluster <- results$clusters
Cite
Please cite our paper if you use this code in your own work:
Song, Q., Su, J., Miller, L.D. and Zhang, W., 2021. scLM: automatic detection of consensus gene clusters across multiple single-cell datasets. Genomics, Proteomics & Bioinformatics, 19(2), pp.330-341.