Merge pull request #87 from openproblems-bio/jalil

Restructuring process data folder + fixing bugs in geneformer

Author: janursa

scripts/run_all.sh

@@ -1,11 +1,11 @@
 set -e
 
-datasets=( 'op' ) #'replogle' 'op' 'nakatake' 'adamson' 'norman'  'xaira_HEK293T' 'xaira_HCT116'  'parsebioscience' 'ibd' '300BCG'
+datasets=( 'replogle') #'replogle' 'op' 'nakatake' 'adamson' 'norman'  'xaira_HEK293T' 'xaira_HCT116'  'parsebioscience' 'ibd' '300BCG'
 
-run_local=false # set to true to run locally, false to run on AWS
+run_local=true # set to true to run locally, false to run on AWS
 
-run_grn_inference=true
-run_grn_evaluation=false
+run_grn_inference=false
+run_grn_evaluation=true
 run_download=false
 
 

scripts/run_grn_evaluation.sh

@@ -65,7 +65,7 @@ if [ -z "${DATASET:-}" ]; then
 fi
 
 num_workers=10
-metric_ids="[regression_2, ws_distance, sem]" #regression_1, regression_2, ws_distance
+metric_ids="[regression_2, ws_distance, sem, tf_recovery, tf_binding, replica_consistency]" #regression_1, regression_2, ws_distance
 RUN_ID="${DATASET}_evaluation"
 models_folder="${DATASET}/"
 apply_skeleton=false
@@ -109,8 +109,6 @@ param_list:
 HERE
 fi
 
-# Write YAML header
-
 append_entry() {
   local grn_name="$1"
   local prediction="$2"
@@ -121,7 +119,7 @@ append_entry() {
       layer_=$layer
   fi
   cat >> "$param_local" << HERE
-  - id: ${reg_type}_${grn_name}_${dataset}
+  - id: ${grn_name}_${dataset}
     metric_ids: ${metric_ids}
     evaluation_data: ${resources_dir}/grn_benchmark/evaluation_data/${dataset}_bulk.h5ad
     tf_all: ${resources_dir}/grn_benchmark/prior/tf_all.csv
@@ -134,12 +132,12 @@ append_entry() {
     layer: $layer_
 
 HERE
-
   # Additional fields for specific datasets
   if [[ "$dataset" =~ ^(norman|replogle|adamson|xaira_HCT116|xaira_HEK293T)$ ]]; then
     cat >> "$param_local" << HERE
     ws_consensus: ${resources_dir}/grn_benchmark/prior/ws_consensus_${dataset}.csv
     ws_distance_background: ${resources_dir}/grn_benchmark/prior/ws_distance_background_${dataset}.csv
+    evaluation_data_de: ${resources_dir}/grn_benchmark/evaluation_data/${dataset}_de.h5ad
 HERE
   fi
 }

scripts/run_grn_inference.sh

@@ -118,8 +118,7 @@ if [[ "$DATASET" =~ ^(replogle|parsescience|xaira_HEK293T|xaira_HCT116)$ ]]; the
   append_entry "$DATASET" "[pearson_corr, negative_control, positive_control, grnboost, ppcor, portia, scenic]"
   append_entry "$DATASET" "[scprint]" "true"
 elif [ "$DATASET" = "op" ] || [ "$DATASET" = "ibd" ]; then
-  # append_entry "$DATASET" "[pearson_corr, negative_control, positive_control, grnboost, ppcor, portia, scenic, scprint, figr, scenicplus, celloracle, granie, scglue]"
-  append_entry "$DATASET" "[ scenicplus, celloracle, geneformer]"
+  append_entry "$DATASET" "[pearson_corr, spearman_corr, negative_control, positive_control, grnboost, ppcor, portia, scenic, scprint, geneformer, scgpt, figr, scenicplus, celloracle, granie, scglue]"
 
 else
   append_entry "$DATASET" "[pearson_corr, negative_control, positive_control, grnboost, ppcor, portia, scenic, scprint]"

scripts/sync_resources.sh

@@ -26,4 +26,7 @@ set -e
 # aws s3 sync  resources/extended_data/ s3://openproblems-data/resources/grn/extended_data --delete
 # aws s3 sync resources/results/experiment  s3://openproblems-data/resources/grn/results/experiment --delete 
 # aws s3 sync resources_test  s3://openproblems-data/resources_test/grn/ --delete 
-aws s3 sync  resources/grn_benchmark/ground_truth s3://openproblems-data/resources/grn/grn_benchmark/ground_truth 
+# aws s3 sync  resources/grn_benchmark/ground_truth s3://openproblems-data/resources/grn/grn_benchmark/ground_truth 
+aws s3 sync  resources/grn_benchmark/evaluation_data s3://openproblems-data/resources/grn/grn_benchmark/evaluation_data --delete
+
+# aws s3 sync   s3://openproblems-data/resources/grn/grn_benchmark/ground_truth resources/grn_benchmark/ground_truth --no-sign-request

src/methods/geneformer/helper.py

@@ -23,7 +23,7 @@
 logger = logging.getLogger(__name__)
 
 
-def tokenize_data(nproc, model_details=None, gene_median=None, token=None, gene_mapping_file=None, tokenized_dir=None):
+def tokenize_data(nproc, temp_dir, model_details=None, gene_median=None, token=None, gene_mapping_file=None, tokenized_dir=None):
     """Tokenize data with required parameters"""
     if not all([model_details, gene_median, token, gene_mapping_file, tokenized_dir]):
         raise ValueError("Missing required parameters for tokenization")
@@ -37,7 +37,7 @@ def tokenize_data(nproc, model_details=None, gene_median=None, token=None, gene_
     )
 
     tokenizer.tokenize_data(
-        "/tmp/geneformer/", tokenized_dir, "tokenized", file_format="h5ad"
+        temp_dir, tokenized_dir, "tokenized", file_format="h5ad"
     )
 # extract embeddings
 def get_embs(
@@ -970,6 +970,7 @@ def tryParallelFunction(func, label, **kwargs):
 
 def compute_geneformer_network(
     adata,
+    temp_dir,
     forward_batch_size=4,
     max_ncells=1000,
     n_processors=20,
@@ -990,19 +991,19 @@ def compute_geneformer_network(
     ]
     adata.obs["n_counts"] = adata.X.sum(1)
     # Create the geneformer folder if it doesn't exist
-    geneformer_folder = "/tmp/geneformer"
+    geneformer_folder = f"{temp_dir}/geneformer"
     if not os.path.exists(geneformer_folder):
         os.makedirs(geneformer_folder)
-    adata.write_h5ad("/tmp/geneformer/to_token.h5ad")
+    adata.write_h5ad(f"{temp_dir}/geneformer/to_token.h5ad")
 
     genelist = [gene_mapping_dict[u] for u in adata.var.index]
 
-    tokenized_data_path = "/tmp/geneformer/tokenized_data.dataset"
+    tokenized_data_path = f"{temp_dir}/geneformer/tokenized_data.dataset"
     if os.path.exists(tokenized_data_path):
         shutil.rmtree(tokenized_data_path)
 
     # Note: This would need proper model_details, gene_median, gene_mapping_file parameters
-    # tryParallelFunction(tokenize_data, "Tokenizing data")
+    tryParallelFunction(tokenize_data, "Tokenizing data", temp_dir=geneformer_folder)
 
     embex = EmbExtractor(
         model_type="Pretrained",  # CellClassifier
@@ -1153,6 +1154,7 @@ def main(par):
         ]
         subadata, net = compute_geneformer_network(
             subadata,
+            temp_dir=par["temp_dir"],
             forward_batch_size=par["batch_size"], 
             n_processors=n_processors,
             max_ncells=par["max_cells"],

src/methods/geneformer/run.sh

@@ -0,0 +1,6 @@
+viash run src/methods/geneformer/config.vsh.yaml -- \
+    --rna resources_test/grn_benchmark/inference_data/op_rna.h5ad \
+    --tf_all resources_test/grn_benchmark/prior/tf_all.csv \
+    --prediction output/prediction.h5ad \
+    --temp_dir output/geneformer
+

src/metrics/replica_consistency/script.py

@@ -33,6 +33,14 @@
   method_id = ad.read_h5ad(par['prediction'], backed='r').uns['method_id']
   dataset_id = ad.read_h5ad(par['evaluation_data'], backed='r').uns['dataset_id']
 
-  output = main(par)
-  
+  try:
+    output = main(par)
+  except Exception as e:
+    print({'error': str(e)})
+
+    output = pd.DataFrame({
+        'key': [None],
+        'value': [None],
+    })
+
   format_save_score(output, method_id, dataset_id, par['score'])

…trics/experimental/sem/config.novsh.yaml src/metrics/sem/config.vsh.yamlsrc/metrics/experimental/sem/config.novsh.yaml renamed to src/metrics/sem/config.vsh.yaml src/metrics/experimental/sem/config.novsh.yaml renamed to src/metrics/sem/config.vsh.yaml

@@ -21,7 +21,7 @@ engines:
     __merge__: /src/api/base_requirements.yaml
     setup:
       - type: python
-        packages: [ tqdm_joblib==0.0.5]
+        packages: [ ]
 runners:
   - type: executable
   - type: nextflow

src/metrics/experimental/sem/helper.py src/metrics/sem/helper.pysrc/metrics/experimental/sem/helper.py renamed to src/metrics/sem/helper.py

@@ -42,7 +42,7 @@ for dataset in "${datasets[@]}"; do
         fi
 
         echo -e "\nProcessing method: $method\n"
-        python src/metrics/experimental/sem/script.py \
+        python src/metrics/sem/script.py \
             --prediction "$prediction" \
             --evaluation_data "$evaluation_data" \
             --score "$score"