search

WWXkenmo / ENIGMA

A fast and accurate deconvolution algorithm based on regularized matrix completion algorithm (ENIGMA)
MIT License
27 stars 6 forks source link

Error: no slot of name "result_CSE_normalized" for this object of class "ENIGMA" #4

Closed Ci-TJ closed 2 years ago

Ci-TJ commented 2 years ago

Hi! Thank you for developing ENIGMA. I tried to reproduce the code of example, but failed. I replaced "result_CSE_normalized" with "result_CSE", but I didn't get the good results as your description. Is the problem with the version of ENIGMA? Rplot01

# loading the packages we need
suppressPackageStartupMessages(library(scater))
suppressPackageStartupMessages(library(splatter))
suppressPackageStartupMessages(library(ENIGMA))
suppressPackageStartupMessages(library(cowplot))

## Simulate CSE
params <- newSplatParams(seed=2022)
params <- setParams(params, update = list(nGenes = 8368, batchCells=200))
sim.groups.c1 <- splatSimulateGroups(params,
                                     group.prob = c(0.25,0.25,0.25,0.25),
                                     de.prob = c(0.3),
                                     verbose = FALSE)
sim.groups.c1 <- logNormCounts(sim.groups.c1)
sim.groups.c1 <- runPCA(sim.groups.c1)
sim.groups.c1 <- runTSNE(sim.groups.c1,dimred="PCA",n_dimred=5)
p_c1 <- plotTSNE(sim.groups.c1, colour_by = "Group",point_size=3)   
label.c1 <- sim.groups.c1$Group

params <- newSplatParams(seed=1002)
params <- setParams(params, update = list(nGenes = 8368, batchCells=200))
sim.groups.c2 <- splatSimulateGroups(params,
                                     group.prob = c(0.5,0.5),
                                     de.prob = c(0.01),
                                     verbose = FALSE)
sim.groups.c2 <- logNormCounts(sim.groups.c2)
sim.groups.c2 <- runPCA(sim.groups.c2)
sim.groups.c2 <- runTSNE(sim.groups.c2,dimred="PCA",n_dimred=5)
p_c2 <- plotTSNE(sim.groups.c2, colour_by = "Group",point_size=3)   
label.c2 <- sim.groups.c2$Group

params <- newSplatParams(seed=1004)
params <- setParams(params, update = list(nGenes = 8368, batchCells=200))
sim.groups.c3 <- splatSimulateGroups(params,
                                     group.prob = c(0.5,0.5),
                                     de.prob = c(0.1),
                                     verbose = FALSE)
sim.groups.c3 <- logNormCounts(sim.groups.c3)
sim.groups.c3 <- runPCA(sim.groups.c3)
sim.groups.c3 <- runTSNE(sim.groups.c3,dimred="PCA",n_dimred=5)
p_c3 <- plotTSNE(sim.groups.c3, colour_by = "Group",point_size=3)   
label.c3 <- sim.groups.c3$Group

p_true <- plot_grid(p_c1,p_c2,p_c3,nrow=1)
p_true <- p_true + ggtitle("Ground Truth") +
    theme(plot.title = element_text(hjust = 0.5))
p_true

################
sim.groups.c1 <- counts(sim.groups.c1) %*% diag(10000/colSums(counts(sim.groups.c1)))
sim.groups.c2 <- counts(sim.groups.c2) %*% diag(10000/colSums(counts(sim.groups.c2)))
sim.groups.c3 <- counts(sim.groups.c3) %*% diag(10000/colSums(counts(sim.groups.c3)))

colnames(sim.groups.c2) <- colnames(sim.groups.c3) <- colnames(sim.groups.c1) <- paste0("Sample-",1:ncol(sim.groups.c1))

mat <- list(sim.groups.c1,sim.groups.c2,sim.groups.c3)
names(mat) <- c("ct1","ct2","ct3")
H1_array <- array(0,
                  dim = c( 3,
                           8368,
                           200))

### construct reference profile
H1 <- array(0,
                  dim = c( 3,
                           8368))
for(i in 1:3){
  H1_array[i,,] <- mat[[i]]
  H1[i,] <- rowMeans(mat[[i]])
}

k <- 3 # number of cell types
ng <- 8368 # number of genes
p <- 200 # number of samples

###Simulate cell type-specific fractions
cc <- matrix(runif(p*k), ncol=k)
cc <- t(scale(t(cc), center=FALSE, scale=rowSums(cc)))
colnames(cc) <- paste('cellType', 1:ncol(cc), sep="")

##evaluate differential expression genes
##calculate 
G <- NULL
for(i in 1:200){
    G <- cbind(G, t(as.matrix(t(as.matrix(cc[i,])) %*% H1_array[,,i])))
}

noise <- t(matrix(rnorm(p*ng)*15, ncol=ng))
noise[noise<0] <- 0
G <- G + noise
k <- apply(G,1,var) > 10^-8
G <- G[k,]
H1 <- H1[,k]
H1_array <- H1_array[,k,] 

rownames(G) <- colnames(H1) <- paste0("gene",c(1:nrow(G)))
colnames(G) <- paste0("Sample",1:ncol(G))
rownames(H1) <- colnames(cc)

raw <- SingleCellExperiment(assays=list(logcounts = G))
raw$Group <- label.c1                
raw <- raw
raw <- runPCA(raw)
p <- plotPCA(raw, colour_by = "Group")
#raw <- runUMAP(raw,dimred="PCA",n_dimred=5)
#p_raw_umap <- plotUMAP(raw, colour_by = "Group",point_size=3)  
raw <- runTSNE(raw,dimred="PCA",n_dimred=5)
p_raw1 <- plotTSNE(raw, colour_by = "Group",point_size=3)

raw$Group <- label.c2
p_raw2 <- plotTSNE(raw, colour_by = "Group",point_size=3)

raw$Group <- label.c3
p_raw3 <- plotTSNE(raw, colour_by = "Group",point_size=3)

p_raw <- plot_grid(p_raw1,p_raw2,p_raw3,nrow=1)
p_raw <- p_raw + ggtitle("Pseudo-Bulk") +
    theme(plot.title = element_text(hjust = 0.5))
p_raw
#############################
get_proportion <- function(X, ref) {
    cat( date(), "Calculating cell type proportion of bulk samples... \n" )
    gene_id = intersect( rownames(X), rownames(ref) )
    X_m = X[gene_id,]
    ref_m = ref[gene_id,]
    ref_m <- apply(ref_m,2,scale)
    theta <- NULL
    coefVec <- NULL
    for(i in 1:ncol(X_m)){
        Exp <- as.matrix(X_m[,i])
        rownames(Exp) <- rownames(ref_m)
        colnames(Exp) <- colnames(X_m)[i]
        Exp <- scale(Exp)

        rlm.o <- rlm(Exp ~ as.matrix(ref_m), maxit = 100)
        coef.v <- summary(rlm.o)$coef[2:(ncol(as.matrix(ref_m)) + 1), 1]
        coefVec <- rbind(coefVec,coef.v)
        coef.v[which(coef.v < 0)] <- 0
        total <- sum(coef.v)
        coef.v <- coef.v/total
        theta <- rbind(theta,coef.v)
    }
    colnames(theta) <- colnames(coefVec) <- colnames(ref_m)
    rownames(theta) <- rownames(coefVec) <- colnames(X_m)
    res <- list()
    res$theta <- theta
    return(res)
}
suppressPackageStartupMessages(library(MASS))
Fra_Simulate <- get_proportion(G, t(H1))

par(mfrow=c(1,3))
plot(Fra_Simulate$theta[,1],cc[,1],main = paste("Cor = ",round(cor(Fra_Simulate$theta[,1],cc[,1]),2),sep=""),xlab="Predict Fractions of CT1",ylab="Ground Truth of CT1")
plot(Fra_Simulate$theta[,2],cc[,2],main = paste("Cor = ",round(cor(Fra_Simulate$theta[,2],cc[,2]),2),sep=""),xlab="Predict Fractions of CT2",ylab="Ground Truth of CT2")
plot(Fra_Simulate$theta[,3],cc[,3],main = paste("Cor = ",round(cor(Fra_Simulate$theta[,3],cc[,3]),2),sep=""),xlab="Predict Fractions of CT3",ylab="Ground Truth of CT3")

suppressPackageStartupMessages(library(TCA))
bulk <- G
time.tca <- system.time({tca.mdl <- tca(X = bulk, W = Fra_Simulate$theta, C1 = NULL, C2 = NULL,
                parallel = TRUE,num_cores=4,max_iters=5)
Z_hat_simulate <- tensor(X = (as.matrix(bulk)), tca.mdl)})

TCA <- SingleCellExperiment(assays=list(logcounts = Z_hat_simulate[[1]]))
TCA$Group <- label.c1                
TCA <- TCA[Fra_Simulate$theta[,1]>0.05]
TCA <- runPCA(TCA)
p <- plotPCA(TCA, colour_by = "Group")
#TCA <- runUMAP(TCA,dimred="PCA",n_dimred=5)
#p_TCA_umap <- plotUMAP(TCA, colour_by = "Group",point_size=3)  
TCA <- runTSNE(TCA,dimred="PCA",n_dimred=5)
p_TCA_tsne1 <- plotTSNE(TCA, colour_by = "Group",point_size=3)  

TCA <- SingleCellExperiment(assays=list(logcounts = Z_hat_simulate[[2]]))
TCA$Group <- label.c2                
TCA <- TCA[Fra_Simulate$theta[,2]>0.05]
TCA <- runPCA(TCA)
p <- plotPCA(TCA, colour_by = "Group")
#TCA <- runUMAP(TCA,dimred="PCA",n_dimred=5)
#p_TCA_umap <- plotUMAP(TCA, colour_by = "Group",point_size=3)  
TCA <- runTSNE(TCA,dimred="PCA",n_dimred=5)
p_TCA_tsne2 <- plotTSNE(TCA, colour_by = "Group",point_size=3)  

TCA <- SingleCellExperiment(assays=list(logcounts = Z_hat_simulate[[3]]))
TCA$Group <- label.c3                
TCA <- TCA[Fra_Simulate$theta[,3]>0.05]
TCA <- runPCA(TCA)
p <- plotPCA(TCA, colour_by = "Group")
#TCA <- runUMAP(TCA,dimred="PCA",n_dimred=5)
#p_TCA_umap <- plotUMAP(TCA, colour_by = "Group",point_size=3)  
TCA <- runTSNE(TCA,dimred="PCA",n_dimred=5)
p_TCA_tsne3 <- plotTSNE(TCA, colour_by = "Group",point_size=3)

p_TCA <- plot_grid(p_TCA_tsne1,p_TCA_tsne2,p_TCA_tsne3,nrow=1)
p_TCA <- p_TCA + ggtitle("TCA") +
    theme(plot.title = element_text(hjust = 0.5))   
p_TCA

############################
G <- as.matrix(G)
H1 <- as.matrix(H1)
egm = create_ENIGMA(bulk = G, ref = t(H1), ref_type = "bulk", meta_ref = as.matrix(colnames(t(H1))))
## Wed Sep 15 11:05:22 2021 Reference from FACS Bulk RNA-seq/microarray.
egm = batch_correct(egm)

## Here we don't need batch effect correction therefore I set the correct bulk expression profile back to the raw!
egm@bulk <- G

## Caculate cell type-specific fractions
egm = get_cell_proportion(egm)

## Using ENIGMA to infer CSE
egm = ENIGMA_L2_max_norm(egm, epsilon=0.001, alpha=0.7,
                        beta=4500,tao_k=0.01,max.iter=1000,verbose=TRUE)
###################################################################
#Error: no slot of name "result_CSE_normalized" for this object of class "ENIGMA"
###################################################################

##########################################
#replaced by following code for the version of ENIGMA
##########################################
#enigma <- egm@result_CSE_normalized[,egm@result_CSE_normalized$cell_type %in% "cellType1"]
enigma <- egm@result_CSE[,egm@result_CSE$cell_type %in% "cellType1"]
enigma$Group <- label.c1                 
enigma <- enigma[Fra_Simulate$theta[,1]>0.05]
enigma <- runPCA(enigma)
p <- plotPCA(enigma, colour_by = "Group")
#enigma <- runUMAP(enigma,dimred="PCA",n_dimred=5)
#p_enigma_umap <- plotUMAP(enigma, colour_by = "Group",point_size=3)    
enigma <- runTSNE(enigma,dimred="PCA",n_dimred=5)
p_enigma_tsne1 <- plotTSNE(enigma, colour_by = "Group",point_size=3)    

#replaced by following code for the version of ENIGMA
#enigma <- egm@result_CSE_normalized[,egm@result_CSE_normalized$cell_type %in% "cellType3"]
enigma <- egm@result_CSE[,egm@result_CSE$cell_type %in% "cellType3"]
enigma$Group <- label.c3                 
enigma <- enigma[Fra_Simulate$theta[,3]>0.05]
enigma <- runPCA(enigma)
p <- plotPCA(enigma, colour_by = "Group")
#enigma <- runUMAP(enigma,dimred="PCA",n_dimred=5)
#p_enigma_umap <- plotUMAP(enigma, colour_by = "Group",point_size=3)    
enigma <- runTSNE(enigma,dimred="PCA",n_dimred=5)
p_enigma_tsne3 <- plotTSNE(enigma, colour_by = "Group",point_size=3)    

#replaced by following code for the version of ENIGMA
#enigma <- egm@result_CSE_normalized[,egm@result_CSE_normalized$cell_type %in% "cellType2"]
enigma <- egm@result_CSE[,egm@result_CSE$cell_type %in% "cellType3"]
enigma$Group <- label.c2                 
enigma <- enigma[Fra_Simulate$theta[,2]>0.05]
enigma <- runPCA(enigma)
p <- plotPCA(enigma, colour_by = "Group")
#enigma <- runUMAP(enigma,dimred="PCA",n_dimred=5)
#p_enigma_umap <- plotUMAP(enigma, colour_by = "Group",point_size=3)    
enigma <- runTSNE(enigma,dimred="PCA",n_dimred=5)
p_enigma_tsne2 <- plotTSNE(enigma, colour_by = "Group",point_size=3)    

p_enigma <- plot_grid(p_enigma_tsne1,p_enigma_tsne2,p_enigma_tsne3,nrow=1)
p_enigma <- p_enigma + ggtitle("ENIGMA") +
    theme(plot.title = element_text(hjust = 0.5))   
p_enigma

sessionInfo()
R version 4.1.0 (2021-05-18)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19042)

Matrix products: default

locale:
[1] LC_COLLATE=Chinese (Simplified)_China.936  LC_CTYPE=Chinese (Simplified)_China.936   
[3] LC_MONETARY=Chinese (Simplified)_China.936 LC_NUMERIC=C                              
[5] LC_TIME=Chinese (Simplified)_China.936    

attached base packages:
[1] parallel  stats4    stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] MASS_7.3-54                 pheatmap_1.0.12             magrittr_2.0.1             
 [4] BiocManager_1.30.16         scater_1.20.1               ggplot2_3.3.5              
 [7] scuttle_1.2.1               TCA_1.2.1                   cowplot_1.1.1              
[10] splatter_1.16.1             SingleCellExperiment_1.14.1 SummarizedExperiment_1.22.0
[13] Biobase_2.52.0              GenomicRanges_1.44.0        GenomeInfoDb_1.28.4        
[16] IRanges_2.26.0              S4Vectors_0.30.2            BiocGenerics_0.38.0        
[19] MatrixGenerics_1.4.3        matrixStats_0.61.0          ENIGMA_0.1.1               

loaded via a namespace (and not attached):
  [1] Rtsne_0.15                ggbeeswarm_0.6.0          colorspace_2.0-2         
  [4] ellipsis_0.3.2            futile.logger_1.4.3       XVector_0.32.0           
  [7] BiocNeighbors_1.10.0      rstudioapi_0.13           farver_2.1.0             
 [10] bit64_4.0.5               AnnotationDbi_1.54.1      fansi_0.5.0              
 [13] splines_4.1.0             sparseMatrixStats_1.4.2   cachem_1.0.6             
 [16] config_0.3.1              nloptr_1.2.2.3            annotate_1.70.0          
 [19] png_0.1-7                 compiler_4.1.0            httr_1.4.2               
 [22] backports_1.4.1           assertthat_0.2.1          Matrix_1.4-0             
 [25] fastmap_1.1.0             limma_3.48.3              BiocSingular_1.8.1       
 [28] formatR_1.11              tools_4.1.0               rsvd_1.0.5               
 [31] gtable_0.3.0              glue_1.6.0                GenomeInfoDbData_1.2.6   
 [34] reshape2_1.4.4            dplyr_1.0.7               gmodels_2.18.1           
 [37] tinytex_0.36              Rcpp_1.0.7                vctrs_0.3.8              
 [40] Biostrings_2.60.2         gdata_2.18.0              nlme_3.1-153             
 [43] DelayedMatrixStats_1.14.3 xfun_0.29                 stringr_1.4.0            
 [46] beachmat_2.8.1            lifecycle_1.0.1           irlba_2.3.5              
 [49] gtools_3.9.2              XML_3.99-0.8              edgeR_3.34.1             
 [52] zlibbioc_1.38.0           scales_1.1.1              RColorBrewer_1.1-2       
 [55] lambda.r_1.2.4            yaml_2.2.1                memoise_2.0.1            
 [58] pbapply_1.5-0             gridExtra_2.3             stringi_1.7.6            
 [61] RSQLite_2.2.9             genefilter_1.74.1         ScaledMatrix_1.0.0       
 [64] checkmate_2.0.0           BiocParallel_1.26.2       rlang_0.4.12             
 [67] pkgconfig_2.0.3           bitops_1.0-7              pracma_2.3.6             
 [70] lattice_0.20-45           purrr_0.3.4               labeling_0.4.2           
 [73] bit_4.0.4                 tidyselect_1.1.1          plyr_1.8.6               
 [76] R6_2.5.1                  generics_0.1.1            nnls_1.4                 
 [79] DelayedArray_0.18.0       DBI_1.1.2                 pillar_1.6.4             
 [82] withr_2.4.3               mgcv_1.8-38               survival_3.2-13          
 [85] KEGGREST_1.32.0           RCurl_1.98-1.5            tibble_3.1.6             
 [88] crayon_1.4.2              futile.options_1.0.1      utf8_1.2.2               
 [91] viridis_0.6.2             locfit_1.5-9.4            grid_4.1.0               
 [94] sva_3.40.0                data.table_1.14.2         blob_1.2.2               
 [97] matrixcalc_1.0-5          digest_0.6.29             xtable_1.8-4             
[100] tidyr_1.1.4               munsell_0.5.0             beeswarm_0.4.0           
[103] viridisLite_0.4.0         vipor_0.4.5
WWXkenmo commented 2 years ago

Thx for using our tool! Yes, its because the version of ENIGMA, in ENIGMA v1.1, we add new slot "result_CSE_normalized" which regress out the effects of cell type proportions to our inferred CSE, this question is because during our model updating, it would latently increase the correlation between each gene and cell type fractions, so if you look at the tSNE plot produced from unnormalized CSE, it would have a "stick" shape, so we need to regard each cell type fraction as some sort of normalize factors. You could install our latest version of ENIGMA in the home page 'ENIGMA_v1.1.tar.gz' through following command

install.packages("ENIGMA_v1.1.tar.gz",repos=NULL, type="source",INSTALL_opts=c("--no-multiarch"))

P.S. The faster and more robust version of ENIGMA is coming soon! hope you could keep attention to our algorithm and give any response or suggestions, which absolutely valuable for us!

WWXkenmo commented 2 years ago

If you use the right version of ENIGMA, the result would like this results I have reproduce this result on a fresh new environment.

Ci-TJ commented 2 years ago

Thanks for your timely reply. I reproduced the code successfully based on your tips. But I encountered some problems with other examples.

  1. ENIGMA_toy2: Error in svd(exp.scale) : infinite or missing values in 'x'

    > egm = ENIGMA_L2_max_norm(egm, epsilon=0.001, alpha=0.8,
+                          beta=10000,tao_k=0.01,max.iter=1000,verbose=TRUE)
Thu Dec 23 19:32:10 2021 Optimizing cell type specific expression profile... 
Ratio ranges from: 4339857.139837 - 4393217.794531
Ratio ranges from: 87613.581326 - 134382.646024
Ratio ranges from: 3160.343952 - 4920.321935
Ratio ranges from: 131.245230 - 212.021252
Ratio ranges from: 7.718480 - 19.369518
Ratio ranges from: 0.376672 - 3.670378
Ratio ranges from: 0.021630 - 0.871688
Ratio ranges from: 0.001492 - 0.216986
Ratio ranges from: 0.000119 - 0.054825
Ratio ranges from: 0.000010 - 0.013953
Ratio ranges from: 0.000001 - 0.003567
Ratio ranges from: 0.000000 - 0.000914
Total loss: 192259347.222029
part1:6121083.91027797 part2:213368509.306481 part3:1544.34558665663
Perform Normalization...Error in svd(exp.scale) : infinite or missing values in 'x'
>

  2. Couldn't find the dataset, deconv.rds. deconv <- readRDS("/mnt/data1/weixu/HiDe/deconv.rds")

  3. Error: 'run.Ted' is not an exported object from 'namespace:TED'

    
> tcga.ted <- TED::run.Ted(ref.dat= H1,
+                          X=t(G)*100,
+                          cell.type.labels= rownames(H1),
+                          tum.key="cellType1",
+                          input.type="GEP",
+                          n.cores=5,
+                          pdf.name="Senescence")
Error: 'run.Ted' is not an exported object from 'namespace:TED'

sessionInfo() R version 4.1.0 (2021-05-18) Platform: x86_64-w64-mingw32/x64 (64-bit) Running under: Windows 10 x64 (build 18363)

Matrix products: default

locale: [1] LC_COLLATE=Chinese (Simplified)_China.936 [2] LC_CTYPE=Chinese (Simplified)_China.936
[3] LC_MONETARY=Chinese (Simplified)_China.936 [4] LC_NUMERIC=C
[5] LC_TIME=Chinese (Simplified)_China.936

attached base packages: [1] parallel stats4 stats graphics grDevices utils datasets [8] methods base

other attached packages: [1] TED_1.1.1 TCA_1.2.1
[3] cowplot_1.1.1 splatter_1.16.1
[5] pheatmap_1.0.12 ENIGMA_0.1.1
[7] scater_1.20.1 ggplot2_3.3.5
[9] scuttle_1.2.1 SingleCellExperiment_1.14.1 [11] SummarizedExperiment_1.22.0 Biobase_2.52.0
[13] GenomicRanges_1.44.0 GenomeInfoDb_1.28.4
[15] IRanges_2.26.0 S4Vectors_0.30.2
[17] BiocGenerics_0.38.0 MatrixGenerics_1.4.3
[19] matrixStats_0.61.0 magrittr_2.0.1

loaded via a namespace (and not attached): [1] backports_1.4.1 spam_2.7-0
[3] plyr_1.8.6 sp_1.4-6
[5] splines_4.1.0 BiocParallel_1.26.2
[7] sva_3.40.0 digest_0.6.29
[9] foreach_1.5.1 matrixcalc_1.0-5
[11] viridis_0.6.2 gdata_2.18.0
[13] fansi_0.5.0 checkmate_2.0.0
[15] memoise_2.0.1 RandomFieldsUtils_0.5.6
[17] ScaledMatrix_1.0.0 config_0.3.1
[19] limma_3.48.3 Biostrings_2.60.2
[21] annotate_1.70.0 gmodels_2.18.1
[23] colorspace_2.0-2 blob_1.2.2
[25] dplyr_1.0.7 RcppArmadillo_0.10.7.5.0 [27] tcltk_4.1.0 crayon_1.4.2
[29] RCurl_1.98-1.5 genefilter_1.74.1
[31] survival_3.2-13 zoo_1.8-9
[33] iterators_1.0.13 glue_1.6.0
[35] gtable_0.3.0 nnls_1.4
[37] zlibbioc_1.38.0 XVector_0.32.0
[39] DelayedArray_0.18.0 BiocSingular_1.8.1
[41] maps_3.4.0 scales_1.1.1
[43] futile.options_1.0.1 DBI_1.1.2
[45] edgeR_3.34.1 Rcpp_1.0.7
[47] RandomFields_3.3.13 viridisLite_0.4.0
[49] xtable_1.8-4 bit_4.0.4
[51] rsvd_1.0.5 dotCall64_1.0-1
[53] animation_2.7 httr_1.4.2
[55] RColorBrewer_1.1-2 ellipsis_0.3.2
[57] pkgconfig_2.0.3 XML_3.99-0.8
[59] farver_2.1.0 locfit_1.5-9.4
[61] utf8_1.2.2 tidyselect_1.1.1
[63] labeling_0.4.2 rlang_0.4.12
[65] reshape2_1.4.4 AnnotationDbi_1.54.1
[67] munsell_0.5.0 tools_4.1.0
[69] cachem_1.0.6 generics_0.1.1
[71] RSQLite_2.2.9 splancs_2.01-42
[73] geoR_1.8-1 stringr_1.4.0
[75] fastmap_1.1.0 yaml_2.2.1
[77] bit64_4.0.5 purrr_0.3.4
[79] KEGGREST_1.32.0 pbapply_1.5-0
[81] nlme_3.1-153 sparseMatrixStats_1.4.2
[83] formatR_1.11 pracma_2.3.6
[85] compiler_4.1.0 rstudioapi_0.13
[87] beeswarm_0.4.0 png_0.1-7
[89] tibble_3.1.6 stringi_1.7.6
[91] futile.logger_1.4.3 fields_13.3
[93] lattice_0.20-45 Matrix_1.4-0
[95] nloptr_1.2.2.3 vctrs_0.3.8
[97] pillar_1.6.4 lifecycle_1.0.1
[99] BiocNeighbors_1.10.0 data.table_1.14.2
[101] bitops_1.0-7 irlba_2.3.5
[103] R6_2.5.1 gridExtra_2.3
[105] vipor_0.4.5 codetools_0.2-18
[107] lambda.r_1.2.4 MASS_7.3-54
[109] gtools_3.9.2 withr_2.4.3
[111] GenomeInfoDbData_1.2.6 mgcv_1.8-38
[113] grid_4.1.0 beachmat_2.8.1
[115] tidyr_1.1.4 DelayedMatrixStats_1.14.3 [117] Rtsne_0.15 ggbeeswarm_0.6.0


By the way, I didn't found R package, scatter. Do you mistyped 'scater' as 'scatter'
Best,
Ci
WWXkenmo commented 2 years ago

let me go straight to your question,

  1. I have fixed the bug, and you could install ENIGMA_v1.2 to see the update
  2. the deconv.rds file is uploaded on our master branch https://github.com/WWXkenmo/ENIGMA/tree/master
  3. this error probably caused by your TED version, mine is TED_1.1, you could read their github to see more details about how to use TED/BayesPrims https://github.com/Danko-Lab/TED !
Ci-TJ commented 2 years ago

@WWXkenmo Hi! For now, it just has the problem with norm.to.one . I couldn't find this function.

exp <- t(norm.to.one(tcga.ted$res$first.gibbs.res$Znkg[,1,]))
Error in h(simpleError(msg, call)) :
  error in evaluating the argument 'x' in selecting a method for function 't': could not find function "norm.to.one"
WWXkenmo commented 2 years ago

Sorry for my late response, norm.to.one is simply the normalize function in TED packages, seem that they have did a lot of changes to their code. Their normalized profile is the "tcga.ted$res$first.gibbs.res$Znkg.norm" I remember, you could directly used their normalized CSE and do the clustering test.

Ci-TJ commented 2 years ago

Sorry for late reply. I changed the codes based on the suggestion from TED developers. Outputs are as below:

exp <- t(tcga.ted$res$first.gibbs.res$Znkg[,1,])
ted <- SingleCellExperiment(assays=list(logcounts = exp))
ted$Group <- label.c4                
ted <- ted
ted <- runPCA(ted)
p <- plotPCA(ted, colour_by = "Group")
#ted <- runUMAP(ted,dimred="PCA",n_dimred=5)
#p_ted_umap <- plotUMAP(ted, colour_by = "Group",point_size=3)  
ted <- runTSNE(ted,dimred="PCA",n_dimred=5)
p_ted_tsne1 <- plotTSNE(ted, colour_by = "Group",point_size=3)

exp <- t(tcga.ted$res$first.gibbs.res$Znkg[,2,])
ted <- SingleCellExperiment(assays=list(logcounts = exp))
ted$Group <- label.c2                
ted <- ted
ted <- runPCA(ted)
p <- plotPCA(ted, colour_by = "Group")
#ted <- runUMAP(ted,dimred="PCA",n_dimred=5)
#p_ted_umap <- plotUMAP(ted, colour_by = "Group",point_size=3)  
ted <- runTSNE(ted,dimred="PCA",n_dimred=5)
p_ted_tsne2 <- plotTSNE(ted, colour_by = "Group",point_size=3)

exp <- t(tcga.ted$res$first.gibbs.res$Znkg[,3,])
ted <- SingleCellExperiment(assays=list(logcounts = exp))
ted$Group <- label.c3                
ted <- ted
ted <- runPCA(ted)
p <- plotPCA(ted, colour_by = "Group")
#ted <- runUMAP(ted,dimred="PCA",n_dimred=5)
#p_ted_umap <- plotUMAP(ted, colour_by = "Group",point_size=3)  
ted <- runTSNE(ted,dimred="PCA",n_dimred=5)
p_ted_tsne3 <- plotTSNE(ted, colour_by = "Group",point_size=3)

overall

I found that my output is basically the same as your result. I will continue to test other examples and wait for a more powerful version of ENIGMA. Best, Qin