update metric values

openproblems-bio · Nov 23, 2024 · df6409b · df6409b
1 parent 543258f
commit df6409b
Showing 1 changed file with 29 additions and 22 deletions.
diff --git a/src/metrics/ks_statistic_sc_features/script.R b/src/metrics/ks_statistic_sc_features/script.R
@@ -9,11 +9,13 @@ requireNamespace("SingleCellExperiment", quietly = TRUE)
 requireNamespace("SummarizedExperiment", quietly = TRUE)
 
 
-packages <- c("SingleCellExperiment", "SummarizedExperiment", "concaveman", "sp", 
-              "Matrix", "methods", "ggplot2", "ggcorrplot", "MuSiC", "fields", 
-              "MCMCpack", "dplyr", "sf", "RANN", "stats", "reshape2", "RColorBrewer", 
-              "scatterpie", "grDevices", "nnls", "pbmcapply", "spatstat", "gtools", 
-              "RcppML", "NMF")
+packages <- c("SingleCellExperiment", "SummarizedExperiment",
+              "concaveman", "sp", "Matrix", "methods",
+              "ggplot2", "ggcorrplot", "MuSiC", "fields",
+              "MCMCpack", "dplyr", "sf", "RANN", "stats",
+              "reshape2", "RColorBrewer", "scatterpie",
+              "grDevices", "nnls", "pbmcapply", "spatstat",
+              "gtools", "RcppML", "NMF")
 
 # Load each package
 lapply(packages, library, character.only = TRUE)
@@ -22,7 +24,8 @@ lapply(packages, library, character.only = TRUE)
 par <- list(
   input_spatial_dataset = "resources_test/spatialsimbench_mobnew/dataset_sp.h5ad",
   input_singlecell_dataset = "resources_test/spatialsimbench_mobnew/dataset_sc.h5ad",
-  input_simulated_dataset = "resources_test/spatialsimbench_mobnew/simulated_dataset.h5ad",  output = "output.h5ad"
+  input_simulated_dataset = "resources_test/spatialsimbench_mobnew/simulated_dataset.h5ad",
+  output = "output.h5ad"
 )
 meta <- list(
   name = "ks_statistic",
@@ -41,7 +44,7 @@ real_prob_matrix <- input_real_sp$obsm[["celltype_proportions"]]
 colnames(real_prob_matrix) <- paste0("ct", seq_len(ncol(real_prob_matrix)))
 
 sim_log_count <- scuttle::normalizeCounts(t(input_simulated_sp$layers[["counts"]]))
-real_log_count <- real_log_count[  , colnames(real_log_count)  %in% rownames(real_prob_matrix) ]
+real_log_count <- real_log_count[ , colnames(real_log_count)  %in% rownames(real_prob_matrix)]
 
 
 sim_sce <- scater::logNormCounts(SingleCellExperiment::SingleCellExperiment(
@@ -55,11 +58,12 @@ sim_log_count <- SummarizedExperiment::assay(sim_sce, "logcounts")
 # build cell type deconvolution in simulated data ## error
 sim_prob_matrix <- CARD_processing(input_real_sp, input_sc)
 
-sim_log_count <- sim_log_count[  , colnames(sim_log_count)  %in% rownames(sim_prob_matrix) ]
+sim_log_count <- sim_log_count[ , colnames(sim_log_count) %in% rownames(sim_prob_matrix)]
 
 feat_types <- c("L_stats", "celltype_interaction", "nn_correlation", "morans_I")
 
-real_scfeatures_result <- scFeatures::scFeatures(real_log_count,
+real_scfeatures_result <- scFeatures::scFeatures(
+  real_log_count,
   sample = rep("sample1", ncol(real_log_count)),
   spatialCoords = input_real_sp$obs[c("row", "col")],
   feature_types = feat_types,
@@ -68,15 +72,18 @@ real_scfeatures_result <- scFeatures::scFeatures(real_log_count,
   spotProbability =  t(real_prob_matrix)
 )
 
-sim_scfeatures_result <- scFeatures::scFeatures(sim_log_count,
-                                                sample = rep("sample1", ncol(sim_log_count)),
-                                                spatialCoords = list( as.numeric( unlist( lapply( strsplit( colnames(sim_log_count) , "x"), `[`, 1) ) )    ,
-                                                                      as.numeric( unlist( lapply( strsplit( colnames(sim_log_count) , "x"), `[`, 2) ) )  )  ,
-                                                feature_types = feat_types,
-                                                type = "spatial_t",
-                                                species = sc_species,
-                                                spotProbability =  t(sim_prob_matrix) )
-
+sim_scfeatures_result <- scFeatures::scFeatures(
+  sim_log_count,
+  sample = rep("sample1", ncol(sim_log_count)),
+  spatialCoords = list(
+    as.numeric(unlist(lapply(strsplit(colnames(sim_log_count) , "x"), `[`, 1))),
+    as.numeric(unlist(lapply(strsplit(colnames(sim_log_count) , "x"), `[`, 2)))
+  ),
+  feature_types = feat_types,
+  type = "spatial_t",
+  species = sc_species,
+  spotProbability =  t(sim_prob_matrix)
+)
 
 ks_statistic_L_stats <- ks::kde.test(x1 = as.numeric(real_scfeatures_result$L_stats), x2 = as.numeric(sim_scfeatures_result$L_stats))
 ks_statistic_celltype_interaction <- ks::kde.test(x1 = as.numeric(real_scfeatures_result$celltype_interaction), x2 = as.numeric(sim_scfeatures_result$celltype_interaction))
@@ -92,10 +99,10 @@ uns_metric_ids <- c(
 )
 
 uns_metric_values <- c(
-  ks_statistic_L_stats,
-  ks_statistic_celltype_interaction,
-  ks_statistic_nn_correlation,
-  ks_statistic_morans_I
+  ks_statistic_L_stats$zstat,
+  ks_statistic_celltype_interaction$zstat,
+  ks_statistic_nn_correlation$zstat,
+  ks_statistic_morans_I$zstat
 )
 
 cat("Writing output AnnData to file\n")