Revision updates

Makefile

@@ -72,6 +72,7 @@ clean:
 	@rm -rf  .nextflow
 	rm -rf tmp
 	find work/ -type f -print0 | xargs -0 -P 8 -n 100 rm -f
+	@rm -rf work
 	@rm -rf plugins
 	@rm -rf .config
 	@rm -f .bash_history

bin/misc_utils/checkers.R

@@ -1,10 +1,54 @@
-# Check if matched sample names, then return those
 check_common_samples <- function(data) {
-  sample_names <- lapply(data$X, rownames) |> unlist() |> unique()
-  matched <- length(sample_names) == length(data$Y)
-  if (matched) {
-    return(sample_names)
+  if (length(data$X) == 0L) {
+    stop("data$X is empty.")
+  }
+
+  sample_names_list <- lapply(data$X, rownames)
+
+  valid_names <- vapply(
+    sample_names_list,
+    function(ids) {
+      !is.null(ids) &&
+        length(ids) > 0L &&
+        !anyNA(ids) &&
+        all(nzchar(ids)) &&
+        anyDuplicated(ids) == 0L
+    },
+    logical(1)
+  )
+
+  if (!all(valid_names)) {
+    stop("Each modality must have non-empty, unique sample row names.")
+  }
+
+  sample_names <- sample_names_list[[1]]
+
+  matched <- vapply(
+    sample_names_list,
+    function(ids) identical(ids, sample_names),
+    logical(1)
+  )
+
+  if (!all(matched)) {
+    stop("Sample names or their order differ between modalities.")
+  }
+
+  if (is.null(data$Y)) {
+    stop("data$Y is missing.")
+  }
+
+  y_count <- if (is.null(dim(data$Y))) {
+    length(data$Y)
   } else {
-    stop("Failed to find match samples, check!")
+    nrow(data$Y)
+  }
+
+  message("Found ", length(sample_names), " samples in each modality.")
+  message("Found ", y_count, " observations in Y.")
+
+  if (length(sample_names) != y_count) {
+    stop("Sample count in X does not match observation count in Y.")
   }
-}
+
+  return(sample_names)
+}

bin/misc_utils/extract_Xy.R

@@ -29,7 +29,7 @@ parseY <- function(Y, verbose = FALSE) {
 }
 
 # Convert MAE to list of X and Y
-extract_Xy <- function(mae, verbose_target=FALSE) {
+extract_Xy <- function(mae, outcome_type="classification", verbose_target=FALSE) {
   # Note need to transpose back to p * n
   # TODO: add check for dimension match
   # COOP LR do not like delayed matrix, so transform it to S3 matrix
@@ -43,11 +43,22 @@ extract_Xy <- function(mae, verbose_target=FALSE) {
   # MAE response would be a dataframe?
   # In simulated data this is always atomic, but in real, they are
   # transformed to dataframe first, so need to pull it
-  y_temp <- mae$response
-  if (is.data.frame(y_temp)) {
-    Y <- parseY(y_temp |> dplyr::pull(response), verbose = verbose_target)
-  } else {
-    Y <- parseY(y_temp, verbose = verbose_target)
+  if (outcome_type == "survival") {
+    #stop("Survival outcome type is not implemented yet")
+    y_df <- SummarizedExperiment::colData(mae) |> as.data.frame()
+    # And get the time and status columns
+    out <- y_df[, c("time", "status")]
+    return(list(X=X, Y=out))
+  } 
+  if (outcome_type == "classification") {
+    y_temp <- mae$response
+    if (is.data.frame(y_temp)) {
+      Y <- parseY(y_temp |> dplyr::pull(response), verbose = verbose_target)
+    } else {
+      Y <- parseY(y_temp, verbose = verbose_target)
+    }
+    return(list(X=X, Y=Y))
   }
-  return(list(X=X, Y=Y))
+  # Otherwise stop
+  stop("Unsupported outcome_type: ", outcome_type)
 }

bin/python_utils/combine_mdata2df.py

@@ -0,0 +1,46 @@
+import mudata
+import pandas as pd
+import numpy as np
+
+
+# def combine_mdata2df(mdata, concat=True):
+
+#     mod_names = list(mdata.mod.keys())
+
+#     X_df = pd.concat( [mdata[k].to_df().add_prefix(f"{k}_") for k in mod_names], axis=1 )
+#     # Also extract the observation df
+#     y_df = mdata[mod_names[0]].obs
+#     # Should contain response column and is of string yes or no
+#     assert y_df["response"].isin(['yes', 'no']).all(), "Column contains values other than 'yes' or 'no'"
+#     # Converting to numeric binary
+#     y_df.loc[:, "response"] = np.where(y_df["response"] == "yes", 1, 0)
+#     # When supplied not to concat meaning y requires other informations more than response
+#     if not concat:
+#         return X_df, y_df
+#     # Merge all DataFrames together, and drop irrevelant meta information
+#     merged_df = pd.concat([ X_df, y_df[["response"]] ] , axis=1)
+#     return merged_df, mod_names
+
+
+def combine_mdata2df(mdata, target_col="response"):
+    # Takes in a mudata and extract X and y components
+    # Get the Xs as a dataframe of combining all modality together columnwise
+    mod_names = list(mdata.mod.keys())
+        # We also add the modality in front of every feature just like "epigenomics_some_feature_name"
+    X_df = pd.concat(
+        [
+            mdata[mod].to_df().add_prefix(f"{mod}_")
+            for mod in mod_names
+        ],
+        axis=1,
+    )
+    # Also extract the observation df
+    y_df = mdata[mod_names[0]].obs.copy()
+
+    assert y_df[target_col].isin(["yes", "no"]).all(), (
+        "Column contains values other than 'yes' or 'no'"
+    )
+    # Converting to numeric binary
+    y_df[target_col] = np.where(y_df[target_col] == "yes", 1, 0)
+
+    return X_df, y_df, mod_names

modules/sklearn/train/resources/usr/bin/load_classifier_class.py → bin/python_utils/load_classifier_class.py

@@ -1,6 +1,11 @@
+# This script is imported for sklearn classifiers usage
+
+
 import scipy.stats as stats
 import importlib
 
+print("This date should be shown as : 2026")
+
 
 def load_classifier_class(model_name, random_state=42, probability=True):
     # Define models with their respective classes, default parameters, and distributions
@@ -9,11 +14,16 @@ def load_classifier_class(model_name, random_state=42, probability=True):
     # https://stackoverflow.com/questions/33843981/under-what-parameters-are-svc-and-linearsvc-in-scikit-learn-equivalent
 
     # Common params
-    C               = stats.loguniform(1e-4, 1e4)
-    min_sample_leaf = stats.randint(1,6)
-    n_estimators    = stats.randint(50, 501)
-    learning_rate   = stats.uniform(0.01, 1.1)
-    max_features    = ["sqrt", "log2", 100, 500, 1000, None]
+    C                   = stats.loguniform(1e-4, 1e4)
+    min_sample_leaf     = stats.randint(1,6)
+    base_n_estimators   = 10
+    n_estimators        = stats.randint(50, 501)
+    learning_rate       = stats.uniform(0.01, 1.1)
+    max_features        = ["sqrt", "log2", 100, 500, 1000, None]
+    max_depth           = 10
+
+
+
     # Dict to store relevant information of sklearn classifiers
 
     # For logistic regression, it has built-in predict proba and coef
@@ -28,70 +38,58 @@ def load_classifier_class(model_name, random_state=42, probability=True):
             "default_params": {"C": 1.0, "kernel": "linear", "random_state": random_state, "probability": probability},
             "params_dist": {"C": C, "kernel": ["linear"] }
     }
-    # Decision Tree classifier, risk of overfitting, hence require pruning of trees
-    decision_tree_dict = {
-            "class_path": "sklearn.tree.DecisionTreeClassifier",
-            "default_params": {"max_depth": 10, "random_state": random_state},
-            "params_dist": {
-                "criterion": ['gini', 'entropy', 'log_loss'],
-                "max_depth": stats.randint(5, 41),
-                "min_samples_leaf": min_sample_leaf,
-                "max_leaf_nodes": [10, 100, 1000, None]
-            }
-    }
     # Random Forest classifier, should in general work better than single decision tree
     random_forest_dict = {
             "class_path": "sklearn.ensemble.RandomForestClassifier",
-            "default_params": {"n_estimators": 10, "max_features": "sqrt", "max_depth": 10, "random_state": random_state, "n_jobs": -1},
+            "default_params": {
+                "n_estimators": base_n_estimators, "max_features": "sqrt", 
+                "max_depth": max_depth, 
+                "random_state": random_state, "n_jobs": -1
+            },
             "params_dist": {
                 "max_features": max_features,
                 "max_leaf_nodes": [10, 100, 1000, None],
                 "min_samples_leaf": min_sample_leaf
             }
     }
 
-    # Gradient Boost classifier, should tune large number of estimators with slow learning rate
-    gradient_boost_dict = {
-            "class_path": "sklearn.ensemble.GradientBoostingClassifier",
-            "default_params": {"loss":'log_loss', "learning_rate":0.1, "n_estimators":10},
-            "params_dist": { 
-                "n_estimators": n_estimators,
-                "learning_rate": learning_rate,
-                "min_samples_leaf": min_sample_leaf,
-                "max_features": max_features
-            }
-    }
     # MLPClassifier is a neural network classifier
     mlp_dict = {
-    "class_path": "sklearn.neural_network.MLPClassifier",
-    "default_params": {"hidden_layer_sizes": (64,), "max_iter": 500,
-                       "random_state": random_state},
-    "params_dist": {"hidden_layer_sizes": [(32,), (64,), (128,), (64, 32)],
-                    "alpha": stats.loguniform(1e-5, 1e-1)}
+            "class_path": "sklearn.neural_network.MLPClassifier",
+            "default_params": {
+                "hidden_layer_sizes": (64,), "max_iter": 500,
+                "random_state": random_state
+            },
+            "params_dist": {
+                "hidden_layer_sizes": [(32,), (64,), (128,), (64, 32)],
+                "alpha": stats.loguniform(1e-5, 1e-1)
+            }
     }
-    # Mimic xgboost
-    hist_gradient_boost_dict = {
-        "class_path": "sklearn.ensemble.HistGradientBoostingClassifier",
-        "default_params": {"random_state": random_state},
-        "params_dist": {"learning_rate": learning_rate,
-                    "max_leaf_nodes": [15, 31, 63],
-                    "min_samples_leaf": stats.randint(5, 30)}
+
+    # Real xgboost
+    xgboost_dict = {
+            "class_path": "xgboost.XGBClassifier",
+            "default_params": {
+                "random_state": random_state ,
+                "n_estimators": base_n_estimators,
+                "objective": "binary:logistic"},
+            "params_dist": { 
+                "learning_rate": learning_rate,
+                "max_depth": stats.randint(2, 9)
+            }
     }
 
     # ========================
     # Lastly merge all together into a single dictionary
     model_info = {
         "Logit": logit_dict,
         "Linear_SVM": linear_svm_dict,
-        # Decision Tree have risk of overfitting, hence require pruning of trees
-        "Decision_Tree": decision_tree_dict ,
-        # RandomForest shuold in general work better than single decision tree
+        # RandomForest should in general work better than single decision tree
         "Random_Forest": random_forest_dict,
-        # GradientBoost should tune large number of estimators with slow learning rate
-        "Gradient_Boost": gradient_boost_dict,
         # MLPClassifier is a neural network classifier
         "MLP": mlp_dict,
-        "Hist_Gradient_Boost": hist_gradient_boost_dict
+        # XGBoost goes here
+        "XGBoost": xgboost_dict
     }
 
     # Check if valid name of model was input

bin/rhelpers.R

@@ -22,3 +22,71 @@ opt2num <- function(opt_chr) {
                                           as.numeric(x), x))
   return(opt)
 }
+
+
+make_parameter_record <- function(value, treatment) {
+  allowed_treatments <- c(
+    "tuned",
+    "fixed",
+    "default",
+    "data-derived"
+  )
+
+  if (!treatment %in% allowed_treatments) {
+    stop(
+      "Unknown parameter treatment: ",
+      treatment
+    )
+  }
+
+  list(
+    value = value,
+    treatment = treatment
+  )
+}
+
+
+write_selected_hyperparameters <- function(
+    dataset_name,
+    method_name,
+    parameters,
+    selection = NULL,
+    output_path = NULL
+) {
+  if (!requireNamespace("jsonlite", quietly = TRUE)) {
+    stop("The jsonlite package is required")
+  }
+
+  if (is.null(output_path)) {
+    safe_method_name <- method_name |>
+      tolower() |>
+      gsub("[^a-z0-9_-]+", "_", x = _)
+
+    output_path <- paste0(
+      safe_method_name,
+      "-",
+      dataset_name,
+      "_selected_hyperparameters.json"
+    )
+  }
+
+  result <- list(
+    dataset = dataset_name,
+    method = method_name,
+    analysis_stage = "model_selection",
+    selection = selection,
+    parameters = parameters
+  )
+
+  jsonlite::write_json(
+    result,
+    path = output_path,
+    pretty = TRUE,
+    auto_unbox = TRUE,
+    null = "null",
+    na = "null",
+    digits = NA
+  )
+
+  invisible(output_path)
+}