Search code examples
rmatrixigraphcorrelationboot

Make a bootstrap for matrix in R with 100 replications (bootstrap correlation matrices)


Here I show the dataset.

library(igraph)
vetor <- c(1, 5, 3, 8, 2, 9, 3, 2:15, 1, 5, 3, 6, 5, 9, 3 )
matrix(vetor, 7, 7)
matrix<-matrix(vetor, 7, 7)
matcor <- cor(matrix, method = "spearman")

matrixnetwork = graph.adjacency(matcor, mode="undirected", weighted = 
TRUE, add.colnames=NULL, diag=FALSE)
plot(matrixnetwork)

I need calculated the bootstrap correlation matrices by randomly permuting the order of values within each Region of Interest, and then computed the values of graphs for each of the bootstrap matrices. I will run subsequent comparison analyzes (t test)

Corr_fun = function(matrix, indices) {
return(cor(matrix, method= "spearman")) }
BIboot= boot(matrix, statistic= Corr_fun, R = 100)

AIBOOT<-BIboot[["t"]]

this form return a matrix 49x100, when I run the the command below, say: Error in graph.adjacency.dense(adjmatrix, mode = mode, weighted = weighted, : not a square matrix.

g_boot <- graph_from_adjacency_matrix(AIBOOT, mode = "undirected", weighted = TRUE)

Furthermore, it returns me with repeated values

enter image description here

I try too:

library(bootstrap)
theta = function(matrix, indices) {
return(as.matrix(matrix)) }  #funcao para centralidade
BOOTSAVE<-bootstrap(matrix, nboot = 100, theta)

here is result: enter image description here

I try too permute for graph, here "matrixnetwork" is from comand for adjacency_matrix, but this form I need do 1 by 1, I need it 100x in a way that makes it easy to measure the graph and make subsequent comparisons:

matrixnetwork2<-permute(matrixnetwork, sample(vcount(matrixnetwork)))

I try this form, seems to work, here "aibi_imuno_BI_1" is my data. I donta put "cm <- cormat[indices, indices]", because this was repeating the same area 2x

N <- 7L set.seed(2023) 
R <- 100L 
BIboot <- vector("list", length = R) 
for(i in seq.int(R)) { indices <- sample(N, replace = TRUE) 
BIboot[[i]] <- cor(aibi_imuno_BI_1, method = "spearman") } 
 g_boot_list <- lapply(BIboot, graph_from_adjacency_matrix, mode = 
"undirected", weighted = TRUE, add.colnames=NULL, diag=FALSE) plot 
 (g_boot_list[[1]])

after correction is this form:

# BOOT 1
N <- 7L
set.seed(2023)
R <- 100L
BIboot <- vector("list", length = R)
for(i in seq.int(R)) {
indices <- sample(N, replace = TRUE)
BIboot[[i]] <- cor(aibi_imuno_BI_1[indices, ], method = "spearman")
}
g_boot_list <- lapply(BIboot, graph_from_adjacency_matrix, mode = 
"undirected", weighted = TRUE, add.colnames=NULL, diag=FALSE)
plot (g_boot_list[[1]])

Here I have 100 bootstrap correlation matrices by randomly permuting

Finally I try:

 #############BOTPERMUTATION FORMA 2 ####################
N <- 7L
set.seed(2023)
R <- 100L
BIPERMUT<- vector("list", length = R)
for(i in seq.int(R)) {
  indices <- sample(N, replace = TRUE)
  BIPERMUT[[i]] <- permute(matrixnetwork, 
 sample(vcount(matrixnetwork)))  
}

plot (BIPERMUT[[1]])
degree(BIPERMUT[[1]], normalized = FALSE, loops = FALSE)

 #################FORMA 3 #######################
 N <- 7L
set.seed(2023)
R <- 100L
BIboot2 <- vector("list", length = R)
for(i in seq.int(R)) {
  indices <- sample(N, replace = TRUE)
  BIboot2[[i]] <- matrixnetwork
}

plot (BIboot2[[1]])
degree(BIboot2[[1]], normalized = FALSE, loops = FALSE)

Solution

  • The posted bootstrap function is not sub-setting the input correlation matrix, it is repeating the same calculation on the entire matrix R times. And the way it is called is wrong too. The question's code is bootstrapping function cor, not the matrix.

    Here is a solution.
    First create a results list BIboot. Then populate the list by bootstrapping the correlations:

    • sample N rows numbers with replacement;
    • create a matrix of those rows and columns;
    • compute the statistic of interest.

    After the correlations list is created, use it to create a list g_boot_list of the R undirected graphs in the final lapply loop.

    library(igraph)
    #> 
    #> Attaching package: 'igraph'
    #> The following objects are masked from 'package:stats':
    #> 
    #>     decompose, spectrum
    #> The following object is masked from 'package:base':
    #> 
    #>     union
    
    N <- 7L
    fakecors <- runif(choose(N, 2))
    cormat <- matrix(nrow = N, ncol = N)
    cormat[lower.tri(cormat)] <- fakecors
    diag(cormat) <- 0
    cormat[upper.tri(cormat)] <- t(cormat)[upper.tri(t(cormat))]
    
    set.seed(2023)
    R <- 10L
    BIboot <- vector("list", length = R)
    for(i in seq.int(R)) {
      indices <- sample(N, replace = TRUE)
      cm <- cormat[indices, indices]
      BIboot[[i]] <- cor(cm, method = "spearman")
    }
    g_boot_list <- lapply(BIboot, graph_from_adjacency_matrix, mode = "undirected", weighted = TRUE)
    

    Created on 2023-10-05 with reprex v2.0.2