version 0.1.0

DESCRIPTION

@@ -0,0 +1,31 @@
+Package: CBnetworkMA
+Type: Package
+Title: Contrast-Based Bayesian Network Meta Analysis
+Version: 0.1.0
+Authors@R: c(
+  person(given = "Garritt L.", family = "Page", email = "page@stat.byu.edu", role = c("aut", "cre", "cph")),
+	person(given = "Andres F.", family = "Barrientos", email = "abarrientos@fsu.edu", role = c("ctb", "cph")),
+	person(given = "S. McKay", family = "Curtis", email = "s.mckay.curtis@gmail.com", role = c("ctb", "cph")),
+	person(given = "Radford M.", family = "Neal", role = c("ctb", "cph")))
+Maintainer: Garritt L. Page <page@stat.byu.edu>
+Description: A function that facilitates fitting three types of models 
+    for contrast-based Bayesian Network Meta Analysis.  The first model is that which
+    is described in Lu and Ades (2006) <doi:10.1198/016214505000001302>.  The other two 
+    models are based on a Bayesian nonparametric methods that permit ties when comparing 
+    treatment or for a treatment effect to be exactly equal to zero. In addition to the 
+    model fits, the package provides a summary of the interplay  between treatment 
+    effects based on the procedure described in Barrientos, Page, and Lin (2023)
+    <doi:10.48550/arXiv.2207.06561>.
+Depends: R (>= 4.2)
+Suggests: igraph
+License: GPL
+Encoding: UTF-8
+LazyData: true
+NeedsCompilation: yes
+Packaged: 2024-03-05 00:49:40 UTC; gpage
+Author: Garritt L. Page [aut, cre, cph],
+  Andres F. Barrientos [ctb, cph],
+  S. McKay Curtis [ctb, cph],
+  Radford M. Neal [ctb, cph]
+Repository: CRAN
+Date/Publication: 2024-03-07 09:40:14 UTC

MD5

@@ -0,0 +1,17 @@
+23fdfc0dd6a0ce23020f989b65325cf9 *DESCRIPTION
+623754f788af21eae362ce754cc9bdc3 *NAMESPACE
+2caa88a889648abfec73727645eb5a9a *R/Network_postprocess.R
+34c45868d92d7e907bc4a8ae0981092a *R/Rutils.R
+97b491834e941c60936c1392b9785608 *R/clique_extract.R
+6b67c476472a4a89936cf1e83918a191 *R/networkMA_wrapper.R
+ceb5ee461f2bd618f6bfec41d86291f9 *data/smokingCessation.RData
+a66569ff55519bcd81c8e58254561131 *man/clique_extract.Rd
+374f0f09787960d6f1d287fd823821f5 *man/networkMA.Rd
+9336d7ad841bb1c421c9873c9541456d *man/network_graphs.Rd
+d03dff69a4379f44e9d94dc4da5aed7c *man/smokingCessation.Rd
+8dadec9b6781567ce27571653819750b *src/Rutil.c
+9e6fa15a1ef455c036b93ff820c22faa *src/Rutil.h
+de75aea2507cac10d14a75c246267c84 *src/init.c
+eebece8e77090ec713c9462d49f03052 *src/matrix.c
+59361ca1e7c67019f745b86ca5cc5d99 *src/matrix.h
+686d8fab4598b1f7640d165d94688db0 *src/networkMA.c

NAMESPACE

@@ -0,0 +1,10 @@
+useDynLib(CBnetworkMA)
+
+export(networkMA)
+export(network_graphs)
+export(clique_extract)
+
+
+
+importFrom("stats", "dnorm")
+importFrom("utils", "combn")

R/Network_postprocess.R

@@ -0,0 +1,57 @@
+# ordmat: MCMC output
+network_graphs = function(ordmat, gamma=c(0.5, 0.75, 0.9, 0.95, 0.99)){
+
+  # gamma is a vector of probability thresholds so that all pairwise edges
+  # in a directed graph have posterior probability of at least gamma
+  #################################
+  # Network with highest post prob
+  ##################################
+
+
+  # Find network with highest prob
+  ordmatcollapse = sapply(ordmat, function(aux)paste(aux, collapse = "|"))
+  Postordmat = sort(table(ordmatcollapse),decreasing = T)
+  prob = Postordmat[1]/sum(Postordmat)
+  Network = ordmat[[which(names(prob)==ordmatcollapse)[1]]]
+
+  mode_graph <- Network
+  post_prob_mode_graph <- as.vector(prob)
+
+  ##################################
+  # Network with pairwise post prob above threshold (Probpair)
+  ##################################
+  # Find network with highest pairwise post prob
+  Network = ordmat[[1]]-ordmat[[1]]
+  pairProbs = ordmat[[1]]-ordmat[[1]]
+  for(i in 1:(ncol(Network)-1))
+  {
+    for(j in (i+1):ncol(Network))
+    {
+      prob = (table(c(sapply(ordmat,function(x) x[i,j]),c(-1,0,1)))-1)/length(ordmat)
+      Network[i,j] = as.numeric(names(which.max(prob)))
+      pairProbs[i,j] = max(prob)
+      # if(max(prob)<=Probpair)
+      #   Network[i,j] = -1111
+    }
+  }
+
+  # Find Network0 which is the coherent graph in the mcmc closest to Network
+  weightl1 = sapply(1:length(ordmat), function(j1) sum(abs(ordmat[[j1]]-Network)*pairProbs))
+  Network0 = ordmat[[sample(which(weightl1==min(weightl1)),1)]]
+
+
+
+  out <- list()
+  cnt <- 1
+  for(Probpair in gamma){
+    # Remove from Network the edges with pairwise probability less than Probpair
+    Network = ifelse(pairProbs <= Probpair, -1111, Network0)
+    Network = ifelse(lower.tri(Network, diag = T), 0, Network)
+    out[[cnt]] <- Network
+    cnt <- cnt+1
+  }
+  names(out) <- gamma
+  list(mode_graph= mode_graph, post_prob_mode_graph=post_prob_mode_graph,out)
+
+}
+

R/Rutils.R

@@ -0,0 +1,60 @@
+
+# Function to find the longest simple path
+find_longest_simple_path <- function(graph) {
+  longest_path <- NULL
+  max_length <- 0
+
+  for (node in igraph::V(graph)) {
+    paths <- igraph::all_simple_paths(graph, from=node)
+    max_path_length <- max(lengths(paths))
+
+    if (max_path_length > max_length) {
+      longest_path <- paths[[which.max(lengths(paths))]]
+      max_length <- max_path_length
+    }
+  }
+
+  return(longest_path)
+}
+
+# Function to find the NMA path within a clique
+path_NMA <- function(v, graph) {
+  # Extract a subnetwork with the given nodes
+  subgraph <- igraph::induced_subgraph(graph, v)
+  igraph::V(subgraph)$label <- igraph::V(graph)$label[v]
+
+  # Find the longest simple path
+  m <- igraph::get.adjacency(subgraph, sparse = FALSE)
+  if (is.null(colnames(m))) colnames(m) <- rownames(m) <- igraph::V(subgraph)$label
+
+  largest_path <- names(sort(rowSums(m), decreasing = TRUE))
+  size_path <- length(largest_path)
+
+  i <- 1
+  aux <- sum(sapply(igraph::all_simple_paths(subgraph, from=which(colnames(m) %in% largest_path[i]),
+                                             to=which(colnames(m) %in% largest_path[i+1])),
+                    length) == 2) ==
+    sum(sapply(igraph::all_simple_paths(subgraph, from=which(colnames(m) %in% largest_path[i+1]),
+                                        to=which(colnames(m) %in% largest_path[i])),
+               length) == 2)
+
+  if (aux) final_path <- paste(colnames(m)[which(colnames(m) %in% largest_path[i])], "=",
+                               colnames(m)[which(colnames(m) %in% largest_path[i+1])])
+  if (!aux) final_path <-  paste(colnames(m)[which(colnames(m) %in% largest_path[i])], "<",
+                                 colnames(m)[which(colnames(m) %in% largest_path[i+1])])
+
+  if (size_path > 2) {
+    for (i in 2:(size_path-1)) {
+      aux <- sum(sapply(igraph::all_simple_paths(subgraph, from=which(colnames(m) %in% largest_path[i]),
+                                                 to=which(colnames(m) %in% largest_path[i+1])),
+                        length) == 2) ==
+        sum(sapply(igraph::all_simple_paths(subgraph, from=which(colnames(m) %in% largest_path[i+1]),
+                                            to=which(colnames(m) %in% largest_path[i])),
+                   length) == 2)
+
+      if (aux) final_path <- paste(final_path, "=", colnames(m)[which(colnames(m) %in% largest_path[i+1])])
+      if (!aux) final_path <- paste(final_path, "<",  colnames(m)[which(colnames(m) %in% largest_path[i+1])])
+    }
+  }
+  final_path
+}

R/clique_extract.R

@@ -0,0 +1,130 @@
+clique_extract = function(ordmat,
+                              type = "Highest_Post_Prob",
+                              clique_size = NULL,
+                              gamma = 0.95,
+                              plot_graph = FALSE){
+  K <- nrow(ordmat[[1]])
+  if(type == "Highest_Post_Prob"){
+    # Find network with highest prob
+    ordmatcollapse = sapply(ordmat, function(aux)paste(aux, collapse = "|"))
+    Postordmat = sort(table(ordmatcollapse),decreasing = TRUE)
+    prob = Postordmat[1]/sum(Postordmat)
+
+      # Plot network
+    Network = ordmat[[which(names(prob)==ordmatcollapse)[1]]]
+    out = cbind(from=1:ncol(Network),to=1:ncol(Network),color=0)
+    for(i in 1:(ncol(Network)-1))
+    {
+      for(j in (i+1):ncol(Network)){
+        if(Network[i,j]==1)
+          out = rbind(out,c(i,j,2))
+        if(Network[i,j]==-1)
+          out = rbind(out,c(j,i,2))
+        if(Network[i,j]==0)
+        {
+          out = rbind(out,c(i,j,1),c(j,i,1))
+        }
+      }
+    }
+
+    # out <- out[,-3]
+    # out <- out[!(out[,1] == out[,2]),]
+
+    # Blue, one direction, orange equal (both direction)
+    mynet <- igraph::graph_from_data_frame(out,directed = TRUE)
+    graph <- mynet
+    igraph::V(graph)$label <- 1:K
+
+    # Plot network
+    if(plot_graph){
+      plot(graph, layout=igraph::layout.circle)
+    }
+    # Find the NMA path for each large clique
+    if(is.null(clique_size)){
+      v <- suppressWarnings(igraph::largest.cliques(graph))
+    } else {
+      # Find the NMA path for smaller clique
+      v <- suppressWarnings(igraph::cliques(graph, min = 2, max = clique_size))
+
+    }
+  }
+  if(type == "Highest_Pairwise_Post_Prob"){
+    # Find network with highest pairwise post prob
+    Network = ordmat[[1]]-ordmat[[1]]
+    pairProbs = ordmat[[1]]-ordmat[[1]]
+    for(i in 1:(ncol(Network)-1))
+    {
+      for(j in (i+1):ncol(Network))
+      {
+        prob = (table(c(sapply(ordmat,function(x) x[i,j]),c(-1,0,1)))-1)/length(ordmat)
+        Network[i,j] = as.numeric(names(which.max(prob)))
+        pairProbs[i,j] = max(prob)
+        # if(max(prob)<=Probpair)
+        #   Network[i,j] = -1111
+      }
+    }
+
+    # Find Network0 which is the coherent graph in the mcmc closest to Network
+    weightl1 = sapply(1:length(ordmat), function(j1) sum(abs(ordmat[[j1]]-Network)*pairProbs))
+    Network0 = ordmat[[sample(which(weightl1==min(weightl1)),1)]]
+
+    Probpair <- gamma ########## THRESHOLD
+    # Remove from Network the edges with pairwise probability less than Probpair
+    Network = ifelse(pairProbs <= Probpair, -1111, Network0)
+    Network = ifelse(lower.tri(Network, diag = TRUE), 0, Network)
+
+    # Plot network
+    out = cbind(from=1:ncol(Network),to=1:ncol(Network),color=0)
+    for(i in 1:(ncol(Network)-1))
+    {
+      for(j in (i+1):ncol(Network)){
+        if(Network[i,j]==1)
+          out = rbind(out,c(i,j,2))
+        if(Network[i,j]==-1)
+          out = rbind(out,c(j,i,2))
+        if(Network[i,j]==0)
+        {
+          out = rbind(out,c(i,j,1),c(j,i,1))
+        }
+
+      }
+    }
+    out
+
+
+
+    # out <- out[,-3]
+    # out <- out[!(out[,1] == out[,2]),]
+
+
+    mynet <- igraph::graph_from_data_frame(out, directed = TRUE)
+    graph <- mynet
+    igraph::V(graph)$label <- 1:K
+    if(plot_graph){
+      plot(graph, layout=igraph::layout.circle)
+    }
+
+
+    # Find the NMA path for each large clique
+    if(is.null(clique_size)){
+      v <- suppressWarnings(igraph::largest.cliques(graph))
+    } else {
+      # Find the NMA path for smaller clique
+      v <- suppressWarnings(igraph::cliques(graph, min = 2, max = clique_size))
+    }
+
+  }
+  cl_list <- rep(0, length(v))
+  count <- 1
+  for(i in 1:length(v)){
+    for(j in 1:length(v)){
+      if (j != i){
+        if(all(v[[i]] %in% v[[j]])){
+          cl_list[i] <- 1
+          break
+        }
+      }
+    }
+  }
+  suppressWarnings(lapply(v[cl_list==0], path_NMA, graph = graph))
+}