Question about NBS: which nodes are part of significant component?

[Jordi H]

Hi all,

I am using BrainGraph to perform the NBS, but there is something which is unclear to me:

after running the NBS function, I get a result like this:

> pairwise_results$`1 vs 4`$components
$observed
Index: <contrast>
   contrast csize      p.perm   Contrast
      <int> <num>       <num>     <fctr>
1:        1    28 0.000999001 Contrast 1

If I understand correctly, this means that one component has been found that significantly differs between the two group I'm comparing, and this component contains 28 nodes. However, what I cannot find out is which nodes form part of this component. Anyone who can tell me how/where I can see that?

Thanks in advance!

Best regards,

Jordi

[Chris Watson]

You can run `make_brainGraphList` on the NBS object and see the vertices that are part of the component there.

Chris

--
You received this message because you are subscribed to the Google Groups "brainGraph-help" group.
To unsubscribe from this group and stop receiving emails from it, send an email to brainGraph-he...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/brainGraph-help/5dca0c47-b55e-499a-b142-214168720137n%40googlegroups.com.

[Jordi H]

Dear Chris,

Thank you for your quick response, apologies for my delayed one.
I tried to run the make_brainGraphList function on my NBS object, but it seems like I have to specify an atlas, which is problematic, because I ran my analysis on a subset of 30 nodes from the HCP-MMP 1.0 atlas. So if I pass that atlas, I get an error because the number of nodes don't match. I tried to make a custom atlas with the right number of nodes, but then I get a different error:

#create fake atlas
my_atlas <- data.frame(name=paste('Region', 1:30), x.mni=rnorm(30), y.mni=rnorm(30), z.mni=rnorm(30),
  lobe=rep(c(1:15), 2),
  hemi=c(rep('L', 15), rep('R', 15)))
my_atlas2 <- as_atlas(my_atlas)

#run NBS
  NBS_res <- NBS(
    SNConMatR_subset,
    covars_subset,
    contrast_matrix,
    con.type = 't',
    X = design_matrix,
    p.init = 0.001,
    N = 1000
  )

There is nothing that indicates the NBS didn't work, for example, NBS_res$p.mat is filled with values that indicate significance. 
(Printing only last two columns:
[,29] [,30] [1,] 1.401819e-07 7.922362e-07 [2,] 6.424397e-11 7.378210e-08 [3,] 6.607323e-12 3.316333e-07 [4,] 2.509797e-09 7.723459e-08 [5,] 5.476912e-08 9.126639e-07 [6,] 5.767420e-06 5.683223e-06 [7,] 2.434526e-06 2.435551e-05 [8,] 3.357912e-11 6.997353e-07 [9,] 9.197361e-11 4.456985e-07 [10,] 1.257340e-09 4.395931e-06 [11,] 1.052310e-09 6.002073e-08 [12,] 1.101879e-05 1.155734e-04 [13,] 2.904585e-09 1.594341e-06 [14,] 2.895869e-06 6.733578e-06 [15,] 2.350167e-13 3.504669e-09 [16,] 1.541276e-08 3.444758e-07 [17,] 3.105240e-08 3.054232e-06 [18,] 4.996555e-10 2.510014e-08 [19,] 4.877029e-07 9.873686e-08 [20,] 9.265072e-11 1.278160e-08 [21,] 2.269630e-07 3.576512e-06 [22,] 1.805899e-08 1.288585e-06 [23,] 1.166274e-08 1.112404e-07 [24,] 7.696265e-13 3.495189e-08 [25,] 1.848075e-11 2.168258e-10 [26,] 1.198847e-11 4.170803e-11 [27,] 7.196360e-10 3.222085e-08 [28,] 1.573183e-11 2.711770e-09 [29,] 0.000000e+00 6.729248e-11 [30,] 0.000000e+00 0.000000e+00 )

But, when I run

make_brainGraphList(NBS_res,my_atlas2) [1] "Start time: 2025-04-04 15:44:09" | | 0%Error in { : task 1 failed - "invalid first argument"

Am I doing something wrong somewhere? Is there a different way to get the connected component?
Do I understand correctly if I think that p.mat contains p-values only of the components that survived the cluster-based correction? If so I suppose I could manually extract them from there, but if there's a different way that would save me some work.

Also, something I noticed, when I run the NBS using f-tests, the $components shows me the connected components:
> pairwise_results_f$`1 vs 3`$components $observed Index: <contrast> contrast csize p.perm Contrast <int> <num> <num> <fctr> 1: 1 30 0.000999001 Contrast 1

But when I run it with t-tests, this stays empty, even when p.mat is filled with non-zero p-values:

> pairwise_results_t$`1 vs 3`$components $observed Index: <contrast> Empty data.table (0 rows and 4 cols): contrast,csize,p.perm,Contrast

Could this be a bug? Or am I just misunderstanding p.mat here?

(This is part of the reason why it would be more work to extract the results manually from p.mat, because now I have to verify if there are more than 1 connected component for each of my comparisons if I want to use the t-test result, whereas otherwise I could get that from the components)

Thanks in advance,

Jordi

Op vrijdag 28 maart 2025 om 12:21:37 UTC+9 schreef Chris Watson:

[Jordi H]

I will just share my code here that I wrote (aka had ChatGPT write) to extract the components automatically, in case any others are struggling.

Note that my results are stored in a list containing multiple groups being compared against each other individually, so anyone with results stored in a different way may need to adjust that

components_f <- list()

for (name in names(pairwise_results_f)) {
  pmat <- pairwise_results_f[[name]]$p.mat
 
  # Make pmat a matrix
  if (length(dim(pmat)) == 3 && dim(pmat)[3] == 1) {
    pmat <- pmat[, , 1]
  }
 
  pmat[lower.tri(pmat)] <- t(pmat)[lower.tri(pmat)]
 
  n <- nrow(pmat)
  visited_nodes <- rep(FALSE, n)
  components <- list()
 
  while (any(!visited_nodes & rowSums(pmat != 0) > 0)) {
    # Start from a node not yet visited with non-zero connections
    current_component <- c()
    to_visit <- which.max(rowSums((pmat != 0) & !visited_nodes))  # best starting point
    queue <- c(to_visit)
   
    while (length(queue) > 0) {
      node <- queue[1]
      queue <- queue[-1]
     
      if (!visited_nodes[node]) {
        visited_nodes[node] <- TRUE
        current_component <- c(current_component, node)
       
        # Add connected nodes that are not yet visited
        connected <- which(pmat[node, ] != 0 & !visited_nodes)
        queue <- unique(c(queue, connected))
      }
    }
   
    components[[length(components) + 1]] <- sort(current_component)
  }

  components_f[[name]] <- components
}

Op vrijdag 4 april 2025 om 22:59:18 UTC+9 schreef Jordi H: