Add cortical only part

Signed-off-by: Liang Zhang <psychelzh@outlook.com>
psychelzh · Mar 9, 2024 · 4eca151 · 4eca151
1 parent d7fcabf
commit 4eca151
Show file tree

Hide file tree

Showing 3 changed files with 53 additions and 16 deletions.
diff --git a/R/constants.R b/R/constants.R
@@ -93,12 +93,16 @@ params_fmri_special <- dplyr::bind_rows(
   )
 )
 
-params_conmat <- tibble::tibble(
+params_conmat <- tidyr::expand_grid(
   xcpd = c(
     "gsr", # with global signal regression
     "no_gsr" # no global signal regression
   ),
-  atlas = "4S256Parcels"
+  atlas = c(
+    "4S256Parcels",
+    # this is the same as 4S256Parcels but without subcortical regions
+    "Schaefer200Parcels"
+  )
 ) |>
   dplyr::filter(xcpd == "gsr")
 config_fmri <- tidyr::expand_grid(

diff --git a/R/prepare_data_neural.R b/R/prepare_data_neural.R
@@ -4,6 +4,9 @@ prepare_meta_time_series <- function(xcpd, session, task, atlas, run) {
     Sys.getenv("ROOT_BIDS_DERIV"),
     sprintf("xcpd_%s", xcpd)
   )
+  if (atlas == "Schaefer200Parcels") {
+    atlas <- "4S256Parcels"
+  }
   extract_bids_meta(
     fs::path(path_xcpd, "xcp_d"),
     fs::path(path_xcpd, "layout"),
@@ -16,10 +19,10 @@ prepare_meta_time_series <- function(xcpd, session, task, atlas, run) {
   )
 }
 
-prepare_data_fc <- function(meta, ...) {
+prepare_data_fc <- function(meta, ..., cortical_only = FALSE) {
   meta |>
     filter(...) |>
-    summarise(fc = list(calc_fc(path)), .by = subject) |>
+    summarise(fc = list(calc_fc(path, cortical_only)), .by = subject) |>
     mutate(
       user_id = data.camp::users_id_mapping[subject],
       .keep = "unused",
@@ -89,8 +92,17 @@ extract_bids_meta <- function(path_bids, path_bids_db, ...) {
     list_rbind()
 }
 
-calc_fc <- function(path) {
-  lapply(path, data.table::fread, na.strings = "n/a") |>
+calc_fc <- function(path, cortical_only = FALSE) {
+  lapply(
+    path,
+    \(path) {
+      ts <- data.table::fread(path, na.strings = "n/a")
+      if (cortical_only) {
+        ts <- ts[, 1:200] # only the first 200 nodes are cortical
+      }
+      ts
+    }
+  ) |>
     data.table::rbindlist() |>
     fisher_cor() |>
     as.dist()

diff --git a/_scripts/prepare_neural.R b/_scripts/prepare_neural.R
@@ -39,14 +39,15 @@ list(
       ),
       tar_target(
         fc,
-        prepare_data_fc(meta_time_series)
+        prepare_data_fc(meta_time_series, cortical_only = cortical_only)
       )
     ),
     config_fc_calc |>
       dplyr::filter(task != "latent") |>
       dplyr::mutate(
         session = lapply(session, \(x) if (x == "0") character() else x),
-        run = lapply(run, parse_digits)
+        run = lapply(run, parse_digits),
+        cortical_only = atlas == "Schaefer200Parcels"
       )
   ),
   tarchetypes::tar_eval(
@@ -60,14 +61,15 @@ list(
     config_fc_calc |>
       dplyr::filter(task == "latent") |>
       dplyr::mutate(
-        call_list_fc = purrr::map2(
-          xcpd, run,
-          \(cur_xcpd, cur_run) {
+        call_list_fc = purrr::pmap(
+          list(xcpd, atlas, run),
+          \(cur_xcpd, cur_atlas, cur_run) {
             config <- config_fc_calc |>
               dplyr::filter(
                 task != "latent",
                 session != "0",
                 xcpd == cur_xcpd,
+                atlas == cur_atlas,
                 if (cur_run == "full") {
                   grepl("^run\\d$", run)
                 } else {
@@ -94,18 +96,37 @@ list(
     names(which(apply(fd_mean, 1, \(x) all(x <= 0.5) && all(!is.na(x)))))
   ),
   tarchetypes::tar_map(
-    dplyr::distinct(params_conmat, atlas),
-    tarchetypes::tar_file_read(
-      atlas_dseg,
+    dplyr::distinct(params_conmat, atlas) |>
+      dplyr::mutate(
+        atlas_store = ifelse(
+          atlas == "Schaefer200Parcels",
+          "4S256Parcels",
+          atlas
+        )
+      ),
+    names = atlas,
+    tar_target(
+      file_atlas_dseg,
       fs::path(
         Sys.getenv("ROOT_BIDS_DERIV"),
         "xcpd_gsr",
         "xcp_d",
         "atlases",
-        sprintf("atlas-%s", atlas)
+        sprintf("atlas-%s", atlas_store)
       ) |>
         fs::dir_ls(regexp = "tsv"),
-      read = read_tsv(!!.x, show_col_types = FALSE, na = "n/a")
+      format = "file_fast"
+    ),
+    tar_target(
+      atlas_dseg,
+      {
+        result <- read_tsv(file_atlas_dseg, show_col_types = FALSE, na = "n/a")
+        if (atlas == "Schaefer200Parcels") {
+          result[1:200, ]
+        } else {
+          result
+        }
+      }
     )
   )
 )