From dcf0f4d00b72aba7c97c1a2396d2495dd9846204 Mon Sep 17 00:00:00 2001
From: Chunmingl <chunming.liu324@gmail.com>
Date: Thu, 8 Jan 2026 00:07:52 -0500
Subject: [PATCH] update method selection for NA twas results

---
 code/pecotmr_integration/twas_ctwas.ipynb | 27 ++++++++++++++++++++++-
 1 file changed, 26 insertions(+), 1 deletion(-)

diff --git a/code/pecotmr_integration/twas_ctwas.ipynb b/code/pecotmr_integration/twas_ctwas.ipynb
index 482b9a5c..0304d597 100644
--- a/code/pecotmr_integration/twas_ctwas.ipynb
+++ b/code/pecotmr_integration/twas_ctwas.ipynb
@@ -576,6 +576,30 @@
     "    library(readr)\n",
     "\n",
     "\n",
+    "    # function to update method selection\n",
+    "    update_twas_method <- function(twas_df){\n",
+    "        twas_df$method_selected_original <- twas_df$is_selected_method\n",
+    "        twas_df_colnames <- colnames(twas_df)\n",
+    "        twas_df$gene_context <- paste0(twas_df$molecular_id, \"_\", twas_df$context)\n",
+    "        rs_update <- do.call(rbind, lapply(unique(twas_df$gene_context), function(gene_context){\n",
+    "            do.call(rbind, lapply(unique(twas_df$gwas_study), function(study){\n",
+    "                df <- twas_df[twas_df$gene_context==gene_context & twas_df$gwas_study==study,,drop=FALSE]\n",
+    "                df_colnames <- colnames(df)\n",
+    "                all_imputable <- df[df$rsq_cv>=${rsq_cutoff} & df$pval_cv<=${rsq_pval_cutoff} & !is.na(df$twas_z),,drop=FALSE]\n",
+    "                if (nrow(all_imputable)<=1 || all(is.na(all_imputable$twas_z))) {\n",
+    "                    return(df) # gene-context pair does not having any other imputable model available.\n",
+    "                } else {\n",
+    "                    message (\"Updating method selection for \", gene_context, \"-\", study, \".\")\n",
+    "                    selected_method_update <- all_imputable[order(-all_imputable$rsq_cv),]$method[1]\n",
+    "                    df$selected_method_update <- df$method==selected_method_update\n",
+    "                    df$is_selected_method <- df$selected_method_update\n",
+    "                    df <- df[, df_colnames]\n",
+    "                    return(df)\n",
+    "                }\n",
+    "            }))\n",
+    "        }))\n",
+    "        return(rs_update[, twas_df_colnames])\n",
+    "    }\n",
     "\n",
     "    # Load metadata and configuration - let these fail if there are issues\n",
     "    if (${\"TRUE\" if rename_column else \"FALSE\"}) {\n",
@@ -778,7 +802,8 @@
     "        final_results <- do.call(rbind, lapply(twas_results_db, function(x) {\n",
     "            if (!is.null(x$twas_result)) x$twas_result[, c(2,1,14:16,3:13)] else data.frame()\n",
     "        }))\n",
-    "        \n",
+    "        # update methods if twas z-score is NA \n",
+    "        final_results <- update_twas_method(final_results)        \n",
     "        message(paste(\"Writing\", nrow(final_results), \"final TWAS results\"))\n",
     "        fwrite(final_results, ${_output[0]:r}, sep = \"\\t\", compress = \"gzip\")\n",
     "    } else {\n",